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Abstract 

Many combinatorial generating functions can be expressed as combinations of sym- 
metric functions, or extracted as sub-series and specializations from such combina- 
tions. Gessel has outlined a large class of symmetric functions for which the resulting 
generating functions are D-finite. We extend Gessel's work by providing algorithms 
that compute differential equations these generating functions satisfy in the case 
they are given as a scalar product of symmetric functions in Gessel's class. Exam- 
ples of applications to A;-regular graphs and Young tableaux with repeated entries 
are given. Asymptotic estimates are a natural application of our method, which we 
illustrate on the same model of Young tableaux. We also derive a seemingly new 
formula for the Kronecker product of the sum of Schur functions with itself. (This 
article completes the extended abstract published in the proceedings of FPSAC02 
under the title "Effective D-Finite Symmetric Functions".) 



Introduction 



A power series in one variable is called differentiably finite, or simply D-finite, 
when it is solution of a linear differential equation with polynomial coeffi- 
cients. This differential equation turns out to be a convenient data structure 
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for extracting information related to the series and many algorithms oper- 
ate directly on this differential equation. In particular, the class of univariate 
D-finite power series is closed under sum, product, Hadamard product, and 
Borel transform, among other operations, and algorithms computing the cor- 
responding differential equations are known (see for instance [34]). Moreover, 
the coefficient sequence of a univariate D-finite power series satisfies a linear 
recurrence, which makes it possible to compute many terms of the sequence 
efficiently. These closure properties are implemented in computer algebra sys- 
tems [24,31]. Also, the mere knowledge that a series is D-finite gives informa- 
tion concerning its asymptotic behavior. Thus, whether it be for algorithmic 
or theoretical reasons, it is often important to know whether a given series 
is D-finite or not, and it is useful to compute the corresponding differential 
equation when possible. 

D-finiteness extends to power series in several variables: a power series is 
called D-finite when the vector space spanned by the series and its deriva- 
tives is finite-dimensional. Again, this class enjoys many closure properties 
and algorithms are available for computing the systems of linear differential 
equations generating the corresponding operator ideals [4,5]. Algorithmically, 
the key tool is provided by Grobner bases in rings of linear differential op- 
erators and an implementation is available in Chyzak's Mgf un package 1 . An 
additional, very important closure operation on multivariate D-finite power se- 
ries is definite integration. It can be computed by an algorithm called creative 
telescoping, due to Zeilberger [46]. Again, this method takes as input (linear) 
differential operators and outputs differential operators (in fewer variables) 
satisfied by the definite integral. It turns out that the algorithmic realization 
of creative telescoping has several common features with the algorithms we 
introduce here. 

Beyond the multivariate case, Gessel considered the case of infinitely many 
variables and laid the foundations of a theory of D-finiteness for symmetric 
functions [9]. He defines a notion of D-finite symmetric series and obtains sev- 
eral closure properties. The motivation for studying D-finite symmetric series 
is that new closure properties occur and can be exploited to derive the D- 
finiteness of usual multivariate or univariate power series. Thus, the main ap- 
plication of [9] is a proof of the D-finiteness for several combinatorial counting 
functions. This is achieved by describing the counting functions as combina- 
tions of coefficients of D-finite symmetric series, which can then be computed 
by way of a scalar product of symmetric functions. Under certain conditions, 
the scalar product of symmetric functions depending on extra parameters is 
D-finite in those parameters, where D-finiteness is that of (usual) multivariate 
power series. Most of Gessel's proofs are not constructive. In this article, we 



This package is part of the algolib library available at http : //algo . inria.fr/ 
packages/. 



2 



give algorithms that compute the resulting systems of differential equations 
for the scalar product operation. Besides Gessel's work, these algorithms are 
inspired by methods used by Goulden, Jackson, and Reilly in [12,13]. Finally, 
Grobner bases are used to help make these methods into algorithms. One 
outcome is a simplification of the original techniques of [12,13]. 

Considering some enumerative combinatorial problem of a symmetric flavor 
and parameterized by a discrete parameter (denoted by k in the examples 
below), it is often so that the enumeration is solved by first forming a scalar 
product of two symmetric functions in k variables. Moreover, in the exam- 
ples envisioned (the enumeration of /c-regular graphs, of /c-uniform tableaux, 
etc.), this scalar product is the specialization to k variables of a scalar prod- 
uct between two "closed form" symmetric functions in infinitely many vari- 
ables. Both symmetric functions are sufficiently well-behaved that nice "closed 
forms" are obtained under specialization, leading to descriptions in terms of 
linear differential operators that are easy to derive. This nice behavior is well 
exemplified by Eq. (5) and Eq. (8) below and is what delimits the scope of 
our method in applications. 

Additionally, our method extends to other scalar products whose associated 
adjunctions satisfy a certain condition of preservation of degree (see Sec- 
tion 9.1), as well as to the Kronecker product of symmetric functions (see 
Section 9.2). 

A very basic example of application of our method is the enumeration of la- 
beled graphs. A finite graph on n vertices labeled with non-negative integers 
ii, . . . ,i n , of respective valencies v±, . . . , v n , is given as a weight the mono- 
mial . . . xl™. This encoding leads to generating functions that are symmet- 
ric series: the set of all finite simple graphs is enumerated by the product 

g ( x ) = j2 n x iXj = u(i+x i x j ), 

GeG (i,j)£E(G) i<j 

as each edge G E(G) is either in the graph or not. This series is obviously 
invariant under renamings of the Xj's, which motivates the involvement of sym- 
metric function theory in the application. Finite simple graphs whose vertices 
all have valency two are called 2-regular graphs. Such a graph contributes to G 
by a term of the form ■ ■ • . Therefore, extracting the sub-series of G with 
same monomials as in the series expansion of IlieN\{o}(l + x J)i another sym- 
metric series, results in the generating series of 2-regular graphs according to 
the same encoding. By symmetry, monomials based on different sets of indices 
ii, . . . , i n of cardinality n share the same coefficient in this extracted series. In 
this spirit, it will be shown in Section 3 that the number of 2-regular graphs 
on n vertices is given as the coefficient of t n in the series 

G 2 (t) = (exp((pi-p 2 )/2-pj/4),exp(t(pl+p 2 )/2)). 
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Here, the scalar product is a scalar product for symmetric functions, to be 
defined in the next section; it implements the coefficient extraction. The vari- 
ables t, pi, and p 2 can be viewed as standard variables, although pi and p 2 
will be assigned the symmetric function interpretation p 1 = x\ + x 2 + • • • , 
P2 = x\ + x\ + • • • . Our purpose in the present paper is to describe scalar 
products of symmetric functions like G 2 (t) by a linear differential equation. 
By our method, Algorithm 1 below calculates that G 2 (t) satisfies the differ- 
ential equation 

2(1 _ t )G' 2 (t) - t 2 G 2 (t) = 0, 

which is easily solved to recover the classical series G 2 (t) = e-VW/y/T^t. 
More details on this calculation as well as similar examples will be given in 
Section 3. In general, the derived differential equation will not admit of such a 
closed form solution. However it is possible to extract asymptotic information 
on the sequence being enumerated directly from this differential equation. This 
will be exemplified in Section 8. 

This article is organized as follows. After recalling the necessary part of Ges- 
sel's work in Section 1, we start by focusing on the special situation when 
a single argument of the scalar product depends on extra parameters. We 
present an algorithm for computing the differential equations satisfied by the 
scalar product in this case in Section 2. The application to the example of 
/c-regular graphs is detailed in Section 3. Then a special case where the algo- 
rithm can be further refined is described in Section 4. We treat a variant of 
Young tableaux where each element is repeated k times in Section 5. (These 
are in bijection with a generalization of involutions [19].) The general form 
of the main algorithm, when both arguments depend on extra parameters, is 
given in Section 6. Termination and correctness of the main algorithms are 
proved in Section 7. Next, in Section 8 we employ our algorithms to derive 
asymptotic estimates of the enumerating sequences of /c-regular graphs for 
k = 1,2,3,4. Following this approach of experimental mathematics, we state 
a conjecture for general k. A discussion on several extensions and applica- 
tions of the method closes the paper in Section 9, including the calculation of 
a seemingly new formula for the Kronecker product of the sum of all Schur 
functions with itself. 



1 Symmetric D-finite Functions 



In this section, we recall the facts we need about symmetric functions, D-finite 
functions, and symmetric D-finite functions. 
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1.1 Symmetric functions 



We first collect basic definitions, notation, and results of the theory of sym- 
metric functions. We refer to [21,34] for further results. 

Symmetric functions are series in the infinite set of variables x\, x 2 , . . . over a 
field K of characteristic 0, subject to a certain invariance under renumberings 
of the variables. The if-algebra A of symmetric functions is formally defined 
as follows. For each positive integer m, the K-vector space consisting of the 
polynomials of K[x\, . . . ,x m ] that are fixed under any permutation of the 
variables is a graded K- algebra G m , the algebra of symmetric polynomials 
in m variables. Here the grading is with respect to the total degree in the 
m variables and it induces a chain of graded surjective homomorphisms 7r m 
from G m+ i onto G m defined by setting x m+ i to 0. Taking the inverse limit 
(a.k.a. projective limit) of the system ({G m }, {n m }) results in the graded K- 
algebra A of symmetric functions. By restriction of the algebras G m and the 
maps 7r m to homogeneous polynomials in a fixed degree n, the inductive limit 
becomes a vector subspace A n of A. We have the relation A = n>o A n . 

We now recall the definitions of the most frequently used bases of the ring A 
and vector spaces A n . Denote by A = (Ai, . . . , A&) a partition of the integer n. 
This means that n = Ai + • ■ ■ + A& and Ai > • • ■ > A& > 0, which we also 
denote A h n. Alternatively, the power notation A = F 1 • • • k Tk for partitions 
indicates that i occurs rj times in A, for % = 1,2, ... , k. Partitions serve as 
indices for the five principal symmetric function families that we use: 

• the homogeneous symmetric functions h\ = h\ 1 ■ ■ ■ h\ k , for h n defined as 
the sum of all monomials of degree n in x 1 , x 2 , ■ ■ ■ , with possible repetition 
(i.e., with any non- negative exponents), 

• the elementary symmetric functions t\ — e\ l ■ ■ • t\ k , for e n defined as the 
sum of all monomials of degree n in x\, x 2 , ■ ■ ■ , with no possible repetition 
(i.e., with exponents or 1, exclusively), 

• the power symmetric functions p\ = p\ x ■ ■ ■ p\ k , for p n defined as the sum 
of the nth power of all variables, 

• the monomial symmetric functions m A = I] cr (r 1 !r 2 ! . . . ) _1 ^m • • • ^(fc)' 
where a ranges over all permutations of the non-negative integers, 

• the Schur symmetric functions S\, whose intuitive definition is in terms of 
the representations of the permutation group S n , and that can alternatively 
be defined as the limit symmetric function when n tends to infinity of the 
determinant of the n x n-matrix with (i, j)-entry h\._ i+ j. 

When the indices are restricted to all partitions of the same positive inte- 
ger n, any of the five families forms a basis for the vector space of symmetric 
polynomials of degree n in Xi,x 2 , .... On the other hand, any of the three 



5 



families indexed by the integers i 6 N, (pi), {hi), and (e^), is algebraically in- 
dependent over Q and generates the algebra A of symmetric functions over K: 
A = K[pi,p 2 , • • • ] = K[hi, h 2 , ■ ■ .] = K[e±, e 2 , . . .]. In this work, we shall focus 
on the basis (pj), as we shall endow A with a differential structure will regard 
to the variables pi. 

Generating series of symmetric functions live in the larger ring of symmetric 
series, K[t] [\pi,P2, •■■]]■ There, we have the generating series of homogeneous 
and elementary functions: 

Hit) = *£h n t n = exp , E(t) =£eX = exp fa-l)^?) . 

1.2 Scalar product and coefficient extraction 

The ring of symmetric series is endowed with a scalar product defined as a 
bilinear symmetric form such that the bases (h\) and (m x ) are dual to each 
other: 

(m^hn) = 6x tli , (1) 
where 5x,^ is 1 if A = /i and otherwise. 

For a partition in power notation, A = l™ 1 • • • k nk , the normalization constant 

z x ■= l n ' ni \ ■ ■ ■ k nk n k \ 
plays the role of the square of a norm of p\ in the following important formula: 

(pa,P m ) = 5x,pZx. (2) 

The scalar product is a basic tool for coefficient extraction. Indeed, if we write 
F(xi,X2, . . . ) in the form J2x fx m x, then the coefficient of x^ 1 • • • x^ k in F is 
fx = (F,h\), by (1). Moreover, when A = l n , the identity h yields 
a simple way to compute this coefficient when F is written in the basis of 
the p's: 

Theorem 1 (Gessel; Goulden & Jackson) Let 9 be the K-algebra homo- 
morphism from the algebra of symmetric functions over K to the algebra K[[t\] 
of formal power series in t defined by 6{p\) = t, 6{p n ) = for n > 1. Then if 
F is a symmetric function, 

0(F) ^^a' 

n=0 U - 

where a n is the coefficient of 
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Gessel also provides an analogue for this theorem when A = l n 2 m and A = 
l n 3 m [9, Theorems 2-4]. Combinations of other degree patterns quickly become 
arduous to write explicitly. 



1.3 Plethysm 

Plethysm is a way to compose symmetric functions, which in the simplest case, 
amounts to simply scaling the indices on the power sums. This inner law of A, 
denoted u[v] for u, v in A, is, for w = Y^x c\P\, defined by the rules [34] 

PnM = C \PnxX 1 PnxX 2 
X 

(au + /3f)M = cm[iu] + (3v [w], (uv)[w] = u[w]v[w], 

where a, (3 in K. For example, consider that w[p n ] = p n [w], and in particular 
that p n [Pm] = Pnxm- Thus, we see that when we write w G A in the power sum 
basis as w — w(jpi,p2, ■ ■ ■ ,Pk, ■ ■ ■), the scaling effect appears on the indices as 

W[p n ] = W(p lxn ,p 2 xn, ■ ■ -,Pkxn, ■ ■ ■)■ 



1.4 D-finiteness of multivariate series 

Recall that a series F G . . . , x n ]} is D-finite in xi, . . . , x n when the set of 

all partial derivatives and their iterates, d ll+ "' +ln F/dx'i ■ ■ ■ dx l ™, spans a finite- 
dimensional vector space over the field K(x±, . . . , x n ). A D-finite description of 
a series F is a set of differential equations whose solutions in any K(xi, . . . , x n )- 
vector space share this property. A typical example of such a set is a system 
of n differential equations of the form 

df d k f 
qi(x)f(x) + q 2 {x) — {x) + ■■■ + qk(x)-^(x) = 0, 

where % ranges over 1, . . . , n, each qj is in K(x±, . . . , x n ) for 1 < j ' < k, and k 
and qj depend on %. Observe that by a theorem of Stafford [2, Chapter 5], any 
D-finite series F admits a D-finite description consisting of only two differential 
equations. However, we do not know how to benefit from this theoretical result 
in our computational setting, and it will be more efficient to compute in a 
systematic way with non-minimal sets. 

The properties of D-finite series we need here are summarized in the following 
theorem. 
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Theorem 2 (1) The set of D-finite power series forms a K -subalgebra of 
\] for the usual product of series; 

(2) If F is D-finite in xi, . . . ,x n then for any subset of variables x ix , . . . ,x ik 
the specialization of F at x^ — • • • — x ik — is D-finite in the remaining 
variables; 

(3) If P is a polynomial in x±, . . . , x n , then exp P(x) is D-finite in x±, . . . , x n ; 

(4) If F and G are D-finite in the variables x\, . . . , x m+n , then the Hadamard 
product F G with respect to the variables xi,...,x n is D-finite in 

X\ , . . . , X m j rn . 

(Recall that the Hadamard product of two series I] a£ Nfe a a u a J2peN k bpu 13 is 
E ae N fe a a b a u a , where u a = u" 1 ■ ■ ■ u^ k .) 

These properties are classical [34] . The first three are elementary, the last one 
relies on more delicate properties of dimension and is due to Lipshitz [20]. 

We note at this point that it is usually simple in applications to provide a 
D-finite description for a D-finite function, as the latter is most often given 
as a polynomial expression in "atomic" D-finite functions, usually well-known 
special functions. Given a table of atomic D-finite descriptions, one bases 
on the closure properties of Theorem 2 above and uses algorithms described 
in [5] in order to derive a D-finite description for the whole expression. In 
our examples, doing this will be straightforward since our functions will be 
exponentials of polynomials. 



1.5 D-finite symmetric functions 



The definition of D-finiteness for series in an infinite number of variables is 
achieved by generalizing Property (2) in Theorem 2: F G i^[[^i, X2, ■ ■ ■ ]] is 
called D-finite in the infinitely many variables Xi if, for any choice of a finite 
set S of positive integers, the specialization to of each for i not in S results 
in a power series that is D-finite, in the classical sense, in the variables Xi 
for % in S. In this case, all the properties in Theorem 2 hold in the infinite 
multivariate case. 

The definition is then tailored to symmetric series by considering the algebra 
of symmetric series as generated over K by the set {pi , p2 , ■ ■ ■ } : a symmetric 
series is called D-finite when it is D-finite in the p^s. 

Property (4) in Theorem 2 has the following very important consequence: 

Theorem 3 (Gessel) Let f and g be elements of K[[t\, . . . ,tk]][[pi,P2, •••]]> 
D-finite in the pi 's and tj 's, and suppose that g involves only finitely many 
of the Pi 's. Then (f, g) is D-finite in the tj 's provided it is well-defined as a 
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power series. 



We return to the example of regular graphs given in the introduction. We 
shall see in Section 3 that the exponential generating series G 2 of 2-regular 
graphs is given as an extraction of coefficients from the generating series G of 
all finite simple graphs in the form G 2 = (G,exp(h 2 t)) and we shall provide 
the explicit representations 



Both G and exp(h 2 t) are clearly D-finite symmetric by the definition above. 
Now, G 2 is equal to the scalar product 



and thus by Theorem 3 the resulting power series is D-finite in t. Note the 
effect of the requirement that g be dependent on finitely many p^s only in 
the theorem — here exp h 2 t depends on p 1 and p 2 only. As a consequence, the 
scalar product extracts those terms from G that are supported by monomials 
in t, pi, and p 2 only. In other words, we can set all pis to in G when i > 2, 
which yields 



This scalar product is between symmetric functions in finitely many pis. 
1.6 Effective D-finite symmetric closures 

Our work consists in making Theorem 3 effective by giving algorithms for 
producing linear differential equations annihilating (f,g). The input to our 
algorithms consists of closed forms for g and the specialization of / in the finite 
number of pi's appearing in g, from which generators of ideals of differential 
operators which annihilate them can then be computed. 

Providing algorithms to manipulate linear differential equations amounts to 
making the closure properties of univariate D-finite series effective; similarly, 
algorithms operating on systems of linear differential operators make the clo- 
sure properties of multivariate D-finite series effective. Our title is thus mo- 
tivated by the fact that our algorithm makes it possible to compute all the 
information on a scalar product that can be predicted from D-finiteness. Note 
that we do not check that the resulting power series is well-defined: our algo- 
rithm merely computes equations that the scalar product series must satisfy 
if it is well-defined. 




and 



Pi +P2 

2 




G 2 (t) = (exp((p?-p 2 )/2-p2/ 4 ),exp(t(p?+p 2 )/2)). 
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In our examples, we make use of symmetric series that are built by plethysm. 
Closure properties are given by Gessel, but our applications require only a sim- 
ple consequence of Property (3) in Theorem 2, namely that if g is a polynomial 
in the p^s, then h[g] and e[g] are D-finite for h = H(l) and e = E(l). 



2 Algorithm for Scalar Product: the Simple Case 

We proceed to give a new algorithm to compute the differential equation satis- 
fied by a scalar product of two D-finite symmetric series under the hypotheses 
of Theorem 3 and with the additional simplifying condition that only one of 
the symmetric series depends on t. When the number of t variables is 1, the 
output is a single differential equation for which existing computer algebra 
algorithms might find a closed-form solution. In most cases however, no such 
solution exists and we are content with a differential equation from which 
useful information can be extracted. 

The basic tool we use here is non-commutative Grobner bases in extensions 
of Weyl algebras. An introduction to this topic can be found in [30]. By W t , 
we denote the Weyl algebra 

W t = K(t 1 ,...,t k ,d tl ,...,d tk ; 

[d tl ,tj] = S id , [ti,tj] = [d ti ,d tj ] = 0, 1 < i, j < k), 

where the bracket [a, b] denotes ab — ba and 5ij is the Kronecker notation. 
This algebra can be identified with the algebra of linear differential opera- 
tors with coefficients that are polynomial in t — t 1; ...,t fc . We correspond- 
ingly denote W p for variables p — p±, . . . ,p n , as well as d t for d tl , . . . , d tk , d p 
for d Pl , . . . , d Pn , etc. For the algorithm, we work in the extension 

W Ptt (t) = K{t) ® K[t] W p , t 

of the Weyl algebra W Pit = W p Cg>x W t in which the coefficients of the differ- 
ential operators are still polynomial in p but rational in t. Suppose F and G 
belong to if [£][[p]] and are D-finite symmetric series as in Theorem 3. In partic- 
ular, they both satisfy systems of linear differential equations with polynomial 
coefficients from K(t)[p\. We can write these equations as elements of W Ptt (t) 
acting on F and G. The set Xp = anmw Ptt (t) F (resp. Xq) of all operators 
of W Pt t(t) annihilating F (resp. G) is then a left ideal of W Pjt (t). Given as 
input Grobner bases of Xp and Xq, our algorithm outputs non-zero elements 
of the annihilating left ideal ann Wt (t) (F,G). 

To combine elements of Xp and X G in a meaningful way we use the adjunction 
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map, denoted o here 2 , defined for an operator P G W p by imposing the 
relation (P ■ F,G) — (F, P° ■ G) for all series F and G. As a consequence, we 
have the relation (PQ)° = Q P° and the adjoint P° is computed formally 
from pi = id Pi and 6£. = p^/i; in particular (pid Pi y = pj<9 Pi [21]. This makes 
the adjunction map an involution as well as an algebra anti-automorphism 
of W p . Note that, although adjunction extends to W p (t) by setting t\ = ti, 
no adjoint for the c\ can be defined in any consistent way. Assume that an 
adjoint d£ existed. For reasons to be explained later, this adjoint has to be of 
the form ad ti + f3ti + 7 for complex constants a, /3, 7, with aj3 7^ 0. Now, for 
any series F and G we have (d ti • F,G) = (F,d%. ■ G^. Choose any non-zero 
series F independent of U\ then by the method of variation of parameters for 
series, one finds a series G satisfying df.-G — F. Upon evaluation, we obtain 
= (F, F) 7^ 0, a contradiction. 

We now proceed to outline the algorithm for the simple case, meaning that 
from this point on we elect to have F G i.e., F independent of t. 

The condition on F that it does not involve t implies that • F — for % 
from 1 to k. We can use this fact to simplify our calculations. In this case, 
we consider a different annihilator, &nnw p F, hereafter denoted J F . Note that 
J F =I F n W p . 

This allows us to determine the action of combinations of P G J F and Q e Xq- 
For example, given any S G W p , T G W p>t (t), and U G W t (t), 

(F, {P°SU + TQ) ■ G) = (S°P -F,U-G) + (F, TQ ■ G) = 0. 

It follows that, if we can find a combination such that J2j PjSjUj + J2j TjQj = 
R G W t , we have = (F, R ■ G) = R ■ (F, G). Note that each P^Sj is an 
element of J F while each TjQj is an element of Xq. Therefore, we conduct our 
search for an element of annyK t (F, G) by determining a non-zero element of 
(j F W t {t) + X G ) f]W t . We shall prove in Section 7.1 that such an element exists. 
Basically, the goal of our algorithms is to compute sufficiently many non-zero 
elements of {j F W t (t) + J G ) fl W t so as to generate a D-finite description of 
the scalar product. 

Note, however, that while X G is a left W Pjt (t) ideal, J F W t (t) is a right W Ptt {t)- 
ideal and the sums P + Q for P G J F W t (t) and Q G X G do not form an 
ideal. This problem is very similar to the problem of creative telescoping: 
given an ideal X C W Ptt (t), the aim in the first step of this method is to 
determine an element of d p W Ptt (t) + X that does not involve p. There also, 
dpW Pt t(t) := J2j 9 Pj W Ptt (t) is a right ideal. The algorithm we present thus bears 
a non- fortuitous resemblance with that of [37]: in this reference, truncations 
of the left ideal X and of the right ideal d p W Pjt (t) at a given total degree in 
p, d p , d t are recombined linearly, this for higher and higher truncation degrees 

2 Macdonald denotes the adjunction operator by _L. 
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until the corresponding truncation of the intersection (d p W Pit (t) + JJ C]W t is 
non-trivial. In our situation, we determine truncations of the left ideal Xq and 
the right ideal J F W t (t) at a given truncation order, recombine those two vector 
spaces linearly, and iterate over higher and higher truncation orders until the 
corresponding truncation of ^J F W t (t) +Xg) H W t is a D-finite description. 

To some extent, the approach of the present paper also shares features with 
that in [27]. However, this reference focuses on determining a bound on a 
truncation order that permits to compute generators of an intersection L = 
(dpW P:t +I^J nW t for a given ideal / of W Pjt , and also generators for a whole free 
resolution of L. From there, the cohomology groups of the module-theoretic 
integral Wt/L of the quotient module W Pt t/I are derived. Roughly speaking, 
we are not concerned here with more than the first cohomology group, and 
furthermore, we treat the similar but different problem for ideals of W Pt t(t) 
and intersections in W t (t). 

Being a module over Wt(t), the sum J F Wt(t) +Xg is a vector space over K(t). 
It is this second structure that is adapted to our method. We could try using 
the module structure in this section, but this would not generalize to the case 
when also F depends on t. The idea is to use K(t)-\mear algebra in the vector 
space structure to eliminate the d Pi and pi. Roughly speaking, we incrementally 
generate lines in a matrix corresponding to generators of J F Wt(t) + Xq, and 
perform Gaussian elimination to remove the monomials involving p and d p . 

The main loop of the algorithm considers monomials of increasing degree with 
respect to any ordering on the monomials in p, d p ,d t . We use the notation ^ 
to denote the monomial comparison associated with this ordering. We reduce 
each monomial a with respect to (the Grobner bases for) X F and Xq. Note that 
the chosen monomial ordering is the same for both Xq and X F . Equivalently, 
the remainder of the reduction of a monomial a with respect to X F can be 
viewed as the adjoint of the remainder of the reduction of a° with respect 
to X F . However, to reflect the fact that adjunction modifies the variables, 
when reducing with respect to Xp we need to use a different order, specifically, 
the ordering < defined by (3i ^ fl 2 on W p if and only if (5\ -< In our 
implementation, we use the ordering DegRevLex(<9 p > p > d t ) which sorts 
by total degree first, breaking ties by a reverse lexicographic order on the 
variables. The order ^ is then DegRevLex(p > d p ). 

Once we have computed these values, we add two rows (and for sufficiently 
large a only one column) in a matrix where we perform Gaussian elimination 
to cancel entries corresponding to monomials involving p and d p . 

We now state the algorithm more formally as Algorithm 1, followed by an 
example in the next section. After this example, we describe the modifications 
necessary to handle specific cases more efficiently, and how to treat the general 
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case. The proofs that these algorithms work and terminate are delayed until 
Section 7. 

Algorithm 1 (Scalar Product) 

Input: Symmetric functions F G K[[p]} and G G both D-finite 

in p,t, given by D-finite descriptions in W p and W Ptt (t), respectively. 
Output: A D-finite description of (F,G) in W t (t). 

(1) Determine a Grobner basis Q G for the left ideal ann^ t (t) G with respect 
to any monomial ordering ^, as well as a Grobner basis Qpo for the 
right ideal ann Wp F° with respect to the monomial ordering induced by ■< 
on W p ; 

(2) B := {}; 

(3) Iterate through each monomial a in p,d p ,d t ; 

(a) Write a = (3^ with (3 G W p and 7 G K[d t \; 

(b) a F := (f5-{(3re^g F .)) i; 

(c) a G := a - (ared±g G ); 

(d) Introduce a F and a G as two new elements into B and reduce so as 
to eliminate p, d p ; 

(e) Compute the dimension of the ideal generated by B n W t {t). If this 
dimension is 0, break and output B n W t (t). 

Notice, if m — 1, as is the case in our examples, there is only one variable t, 
and the dimension condition in (3e) is simplified to: 

// B contains a non-zero element P from W t (t), break and return P. 

Note that Step (1) requires to determine both ideals eamw p , t (t) G and ann^ F, 
not just annvK Pi i( P ,t) G and annvy p ( p ) F. In other words, one generally needs 
to pass from a D-finite description {P{\ generating the ideal ann^ ( p ) F as 
J2i W p (p)Pi to a set {Qi} generating the ideal ann Wp F = W p fl ann Wp ( P ) F as 
J2iW p Qi, and similarly for G. The operation of computing such intersections 
is called Weyl closure, in the terminology of [40,41]. It is a non-obvious task, 
owing to the change of module structure (coefficients in W p {p) are replaced 
with coefficients in W p ). Algorithms are provided in [40,41]. 

Sometimes, the input set {Pi} already constitutes a generating set for the 
Weyl closure. In this case, one can skip Step (1) of the algorithm. This is the 
case in our examples. 

The remainder of the reduction with respect to the Grobner basis Q G is a 
multivariate analogue of the remainder of the Euclidean division. It is such 
that for any a, ag = a — (aredC?) belongs to the ideal generated by Q. 
A similar statement holds for Q F . 

For this description we have assumed that Grobner bases could be computed 
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for both left and right ideals. If they can only be computed on one side, say for 
left ideals, then the operators otp can be obtained as follows: first, determine 
the monomial ordering ^ induced by adjunction on W p viewed as a left 
structure from the ordering ^ on W p viewed as a right structure; then, replace 
the Grobner basis Qf° with the Grobner basis Qp for the left ideal annw p F 
with respect to a f is then computed as (j3 — (/3°red^ Q F )^. This way 
we get Q F o = {QfY- 

We represent the basis B as a matrix, with columns indexed by monomials 
in the pj's, the d Pi 's, and the 9 ti 's. Each row in the matrix corresponds to the 
row vector of the coefficients of some element of B with regard to the indexing 
monomial basis. Introducing an element into the basis consists of adding a new 
row at the bottom of the matrix, performing row reduction (also known as 
Gaussian elimination), and then returning the new matrix as the updated 
basis. In practice, B can be handled (not inefficiently) by a Grobner basis 
computation with respect to a monomial ordering that eliminates the p^s and 
the dp/s, performing calculations in the free X[t]-module with a basis the list 
of indexing monomials. 

Finally, some remembering can be done at Step (3b) to avoid reducing the 
same f3 again and again, for different a's involving the same (3. 



3 Example: A>Regular Graphs 

The enumeration of regular graphs has been treated by a number of authors 
[6,9,13,29]. We present it here because of its expository value and as it is the 
simplest in a family of examples. After expressing the problem as a scalar 
product, we describe in detail how our algorithm treats it. We conclude this 
section with an indication of how the scenario may be generalized. 



3.1 A generating series for graphs as a scalar product 

Recall from the introduction that a generating series for the set of all finite 
simple graphs labeled with integers from N \ {0} is 

G ( x ) = 2 II x i x 3 = IE 1 + x i x i\ 

GeG (i,j)eE(G) i<3 

under the encoding that a graph on n vertices ii, . . . ,i n of respective valen- 
cies Vi , . . . , v n contributes a monomial x^ 1 . . . . We can similarly make a 
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generating function for graphs with multiple edges (multigraphs) by 



for an edge of a graph with multiplicity m contributes a monomial x™x™ 
and any non-negative multiplicity is allowed. 

Clearly both G and M are symmetric functions, and in fact, we have the 
relations G = e[e 2 ] and M = h[e 2 ], as determined by a method that we 
discuss in Section 3.4. Both are easily rewritten in terms of the p^s: 



G = exp (XX-I)^? - P2i)/2i^ and M = exp ^ ( ft 2 + p 2l ) /2i) 



(3) 



In any given term, the degree of Xk gives the valency of vertex k. So, for 
example, the coefficient g n of X\---x n in G, hereafter denoted [x± ■ ■ ■ x n ]G, 
gives the number of 1-regular graphs, or perfect matchings on the complete 
graph on n vertices, and in general the coefficient g$ = [x\ - ■ ■ x*]G, also given 
as [rrikn]G, gives the number of A;-regular graphs on n vertices. By virtue of 
Eq. (1), coefficient extraction amounts to a scalar product, and the generating 
function G k (t) of /c-regular graphs is given by 

t n 

G k (t):=Y / 9 l n ] ^ = (G,H k ), where 

Hkif) '-=^h k n 



= E^ = exp(M- (4) 



Now, since h k = H\\-kP\l z \ (where the sum is over all partitions A of k), the 
exponential generating function H k (t) is also exp^t^AhnPA/^A)! an exponen- 
tial in a finite number of p^s. By Property (3) in Theorem 2, this is D-finite. 
Further, as a result of scalar product property (2), we can rewrite Eq. (4) as 

g.m = («p( E (-Dif + T + S © 

\ \ievcn, i<k Zf 1 iodd, i<k^ L j \ \hk ^ / I 

and now by Theorem 3 this generating function G k (t) is D-finite. 

Note how the closed form for G in (3), in infinitely many variables, and the 
closed form for H k (t) in (4), in terms of the /i's, have led to the scalar prod- 
uct (5) between two closed forms, explicitly written in terms of finitely many pi 
for each k. This reduction is what has made the algorithm applicable. 
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3.2 Effective Computation for k = 2 



To illustrate a typical calculation, we calculate G 2 (t), the generating function 
for 2-regular graphs which, according to Eq. (5), is determined by 

G 2 (t) = (e W ((p 2 1 -p 2 )/2-p 2 2 /A),e W (t(p 2 1 +p 2 )/2)). 

Algorithm 1 calculates that G 2 (t) satisfies the differential equation 

2(l-t)G' 2 (t)-t 2 G 2 (t) = 0, 

which is easily solved to find G 2 (t) = e~^ t< - t+2 ' ) / 'y/1 — t. 

In order to appeal to Algorithm 1, set F = exp((p 2 — p 2 )/2 — p 2 ,/^) and 
G = exp(t(pl +p 2 )/2) and determine the Grobner bases Qf and Qg of their 
annihilating ideals respectively: 

Qf = {p 2 + 2d P2 + l,p 1 -d pi } and Q G = {28^ - t,c^_- tp 1: pl+p 2 - 2d t }, 

where Qp is a Grobner basis with respect to the degree reverse lexicographical 
monomial ordering such that p± > p 2 > d Pl > d P2 and Qg is a Grobner 
basis with respect to the degree reverse lexicographical monomial ordering 
such that d pi > d P2 > p± > p 2 > d t . (Leading monomials with respect to the 
monomial ordering are underlined.) Before proceeding, the set Qp is converted 
by adjunction into a Grobner basis Q F with respect to the degree reverse 
lexicographical monomial ordering such that d Pl > d P2 > pi > p 2 : 

g F = {2d P2 +p 2 + i,d PL -p 1 }. 

(The reader should not get confused by the peculiar situation of this example: 
here, adjunction has not changed the polynomials, except for signs, but this 
is only a coincidence.) 

The initial value of B is the empty set. For the sake of the example, we shall 
iterate on monomials a according to the degree reverse lexicographical order 
such that d t > d P2 > p 2 > d Pl > p±, and perform reductions when inserting 
into the basis according to the elimination order sorting first by the degree 
reverse lexicographical order such that d P2 > p 2 > d pi > p 1 , and breaking ties 
by the degree in d t . 

We now briefly sketch the run of the algorithm until a becomes p\d pi and then 
illustrate the steps of the main loop in more details. 

For a — 1 and a = pi, the algorithm inserts no polynomial into the basis B. 
The next iteration of the loop, for a = d Pl , produces a F — d pi — p 1 , which is 
inserted into B as is, and a G = d pi —tpi, whose insertion puts p\ into B. Next, 
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the case a = P2 inserts no polynomial before, for a = d P2 , chf — 2d P2 + P2 + 1 
gets inserted as is, and the insertion of = 2d P2 — t puts P2 + (t + 1) into B. 
The iteration for a = dt has no effect on B. For a = pf, cxf = is not inserted, 
and Oic = p\ + P2 — 2<9< gets inserted in the form p\ — 2d t — (t + 1). 

At this point, the algorithm is about to treat a = p±d Pl and the value of B is 
B = {^-p 1 ,p 1 ,2d El + p 2 + l,p2 + (t + l),pl-2d t -(t + l)}, (6) 

where we have written elements in the order of introduction into the set. In 
matrix notation, the column vector of elements of B reads: 



(p\\ 

o \ d P2 

P2 

1 d pi 

Pi 
d t 

V 1 / 



Here, we have chosen to keep the rows in the order of creation by the algo- 
rithm and to sort the column according to the monomial order used by the 
elimination step. Observe that in this way, no two rows have their left-most 
non-zero entry on the same column: simply reordering rows would put the 
matrix in row echelon form. 

Then, the algorithm computes 
cxf = a — (a red^ Q F ) = a — (a° red^ Qf)° = Pid Pl — pi + 1 



f° 





1 


-1 














1 








2 1 














1 











I 1 











-2 



and 



a G = a - (a rcd ^ g G ) = p x d pi + tp 2 - 2td t . 



(Note that ap is really (<9 Pl — P\)pi, an element of the right ideal generated 
by Q F .) Next, we update B to include these two values. We insert a F into B 
after one reduction, leading to 

B := BU{ Pl d pi -2d t -t}. 

In matrix notation, this insertion adds a new column to the left of the matrix, 
corresponding to the new monomial Pid pi , and one more row at the bottom 
of the matrix, (100000-2 -t). Then the algorithm inserts a G - Its leading 
monomial p\d Pl is already present in B, leading to an initial reduction to tp 2 + 
2(1 — t)d t + t. One final reduction by t times the pre-last element in Eq. (6) 
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2-regular graphs 



0i 

02 



-t 2 

-2t + 2 




3-regular graphs 



4>o 



t 3(t A + 2t 2 -2) 2 
-3(t 10 + Qt 8 + 3t 6 - 6t 4 - 26t 2 + 8) 
-9t 3 (t 4 + 2t 2 -2) 



01 



4-regular graphs 



-t 4 (t 5 + 2t 4 + 2t 2 + 8t - 4) 2 
-4(t 13 + 4t 12 - 16t 10 - 10t 9 - 36t 8 - 220t 7 - 348t 6 

-48t 5 + 200t 4 - 336t 3 - 240i 2 + 416t - 96) 

m 2 (t - l) 2 (t 5 + 2t 4 + 2t 2 + 8t - 4)(t + 2) 2 



Table 1 

Differential equation faG'l + 4>\G' k + 00^ = satisfied by Gfc(i), fc = 2, 3, 4. 
results in the step 

B := BU {2(1 -t)dt-t 2 }. 

The intersection of this and W t (t) is non-trivial, and the algorithm outputs 
2(1 — t)d t — t 2 . We conclude that G 2 (t) satisfies the differential equation 

2{l -t)G' 2 {t) -t 2 G 2 {t) = 0. 

Table 1 summarizes the results by the same algorithm for k = 2,3,4. These 
match with the results in [13]. 



3.3 Efficient enumeration of k -regular graphs 



An efficient procedure for the enumeration of A;-regular graphs is immediately 
derived from the differential equations for the generating series of /c-regular 
graphs collected in Table 1. Indeed, one simply needs to convert the differential 
equation for Gk(t) into a recurrence relation for its coefficients g$ and to de- 
termine sufficiently many starting values g^\ gf \ . . . Then, one can efficiently 
compute g$ for any n by unrolling the recurrence. 

Implementations are available to help with this approach. For example, the 
Maple package gfun 3 by Salvy and Zimmerman [31] contains commands ded- 

3 This package is part of the algolib library, which is available at http://algo 
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icated to the conversion step and the iterative calculations based on a linear 
recurrence. Computations in the case k = 4 result in a recurrence relation of 
order 15 already published by Read and Wormald [29] and can be found as 
a formula accompanying sequence number A005815 in Sloane's encyclopedia 
of integer sequences [32] . From this recurrence relation and initial terms, it is 
then a matter of seconds to compute the exact integer values for hundreds of 
terms in the sequence. 

It should be stressed that this method proves much more efficient than the 
direct computation of the scalar product based on a termwise expansion and 
application of formula (2). For example, Stembridge's implementation in the 
package SF for symmetric function manipulation in Maple [35] already requires 
several minutes to compute the gffl for n up to 15, and becomes unsuitable 
to handle the symmetric functions that would be necessary to obtain g^. 
Far from showing any weakness of SF's general approach, this illustrates the 
computational progress provided by our techniques in the specific setting of 
differentiably finite series. 



3.4 Generalization 

The series given by Eq. (3) is determined combinatorially in a direct fashion 
using the theory of species [1]. This can be extended naturally to handle 
a wider family of combinatorial structures, such as hypergraphs, set covers, 
latin rectangles. For an in-depth treatment, consult [26]. 



4 Hammond Series 

In the example above, it turned out that except for monomials of degree 1, 
we needed only examine the two monomials p\ and p\d pi in order to reach the 
solution. However, depending on the monomial ordering, the algorithm might 
well consider many monomials before it adds the ones that eliminate the pi's 
and <9 Pi 's. The problem becomes far more serious as the number of variables 
and the degree of the monomials increase. It turns out that in the common case 
when the scalar product is of the type (F, Hk(t)) it is possible to modify the 
approach and eliminate the Pi and the d Pi in a more efficient manner using the 
Hammond series 4 (or H-series) introduced by Goulden, Jackson, and Reilly 



. inria.fr/packages/. 

4 In [12, Sec. 3.5] this is referred to as the Gamma series of F. 
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in [13]: for F £ K[\pi,p 2 , ■■■]}, the Hammond series of F is defined as 



where the sum is over all partitions, and if A = l mi • • • k mk then t = t™ 1 ■ ■ ■ t™ k 
and m(X)\ = mi! m 2 ! • • -m*,!. These are very closely related to the Hammond 
operators, defined by Hammond [15] and used extensively by MacMahon [22]. 
A Hammond operator can be described as h\, and thus the Hammond series 
of F with all of the t variables set to 1 results essentially in a sum of Hammond 
operators acting on F. 

Observe that the generating function for /c-regular graphs is 



where the t occurs in position k. This is true for any generating function which 
takes the form (F, Hk(t)} for some F. 

A theorem from [13] is specially useful: Goulden, Jackson, and Reilly's H- 
series theorem states that 7i.{d Pn ■ F) and 7i{p n F) can be expressed in terms 
of the dt t -7i(F)'s. In terms of Grobner bases, this corresponds to introducing 
the additional variables t±, tk (instead of t — tk alone) and work with 
the series Tik(ti, ■ ■ ■ ,tk) = J2\ h\t x /m{\)\ with sum over partitions A whose 
largest part is k (instead of working with the univariate H k (t)). The H-series 
theorem therefore implies that for an appropriate monomial order, there is a 
Grobner basis of the ideal In k of all operators of W Ptt annihilating TCk, with 
elements of the form 



where all the Pi and Qi are polynomials in t, dt- 
The algorithm in this case is as follows. 

Algorithm 2 (Hammond Series) 

Input: An integer k, and F e K[\pi, ■ ■ ■ ,p n ]]- 
Output: A differential equation satisfied by 



where tk is in position k. 

(1) Compute Qf, a Grobner basis for the left ideal Jf annihilating F in W p ; 

(2) Compute Gn k > a Grobner basis of the form (7); 




G fc (f)=W(G)(0,...,0,f,0,...) 



Pi - Pi(t, d t ), d Pi - Qi(t, d t ), 



i — 1, . . . , k, 



(7) 
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(3) For each U G Qf, compute ru G W t as the reduction ofU° by for an 
order which eliminates p, d p . Let Rq be the set of ru 's; 

(4) Fori from 1 to k — 1 eliminate d ti from Ri-\ and set ti — in the resulting 
polynomials; call Ri the new set; 

(5) Return Rk-i- 

As with Algorithm 1, the first step is to determine an annihilating ideal in W p . 
Again, one can possibly first determine a D-finite description and use Weyl 
closure [40,41] to obtain the annihilating ideal. 

After Step (3), all the pi's and <9 Pi 's have been eliminated and Rq contains a 
set of generators of a D-finite W t (t)-ideal annihilating (F,Tl k ). Then, in order 
to obtain differential equations satisfied by the specialization at t\ — ■ ■ ■ — 
tk-i = 0, Step (4) proceeds in order by eliminating differentiation with respect 
to ti and then setting ti = in the remaining operators. 

Note that the Grobner basis of Step (2) can be precomputed for the re- 
quired fc's (although most of the time is actually spent in Step (4)). 

In order to compute the elimination in Step (4), one should not compute a 
Grobner basis for an elimination order, since this would in particular perform 
the unnecessary computation of a Grobner basis of the eliminated ideal. In- 
stead, one can modify the main loop in the Grobner basis computation so that 
it stops as soon as sufficient elimination has been performed or revert to skew 
elimination by the non- commutative version of the extended Euclidean algo- 
rithm as described in [5]. This is the method we have adopted in the example 
session given in Appendix B 5 . 

This calculation is comparatively rapid since the size of the basis is greatly 
reduced. Further, the basis grows smaller as the algorithm progresses, on ac- 
count of setting variables to 0. We can compute the case of 4-regular graphs 
in a second, instead of a couple of minutes using the general algorithm. The 5- 
regular expression requires significantly more computation time, and we could 
not compute it. 

A mathematically equivalent but slightly faster way of performing Step (3) is 
to compute r v by simply replacing each monomial p® 1 ■ ■ -Pn"dpl • • • d^" in U 
with the product Q^P* 1 • • ■ P" 1 - 

In order to explain the relative speed of Algorithm 2, compared to Algorithm 1, 
it needs to be said that the Hammond series approach searches a smaller 
space, which can well result in a differential equation of order higher than 
that obtained by Algorithm I. This occurs, for instance, in the case of 4- 

5 An implementation of the algorithms presented here is available in the Maple 
package ScalarProduct available at http://algo.inria.uqam.fr/ mishna. 
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regular graphs: Algorithm 2 returns a differential equation of order 3 only 
when that returned by Algorithm 1 is of order 2. 



In the same vein, note that the order in which the eliminations are done in 
Step (4) could be changed, possibly leading to a different (but correct) output. 

4-1 Proof of Termination and Correctness 

Termination of Algorithm 2 is obvious. On the other hand, the full proof of 
correctness requires technical results to be proved in Section 7. The following 
corollary articulates a property of D-finite functions in the simple language of 
symmetric functions and D-finite descriptions, and is a corollary of Proposi- 
tion 9 that will be proved independently. 

Corollary 4 Let F and G be D-finite symmetric series in X[[pi, . . . ,p n }} an d 
K[ti, . . . ,tk][[pi, . . . ,p n }}, respectively, with corresponding annihilators Jp C 
W p andXc C W Ptt (p,t). Under these conditions, the vector space 



is a D-finite symmetric series by Theorem 2 since each hj is a finite combination 



of pi, p n . Thus, / = H(F)(t 1 , ... ,t k ) = (H k (ti, ■ ■ -,t k ),F) is a D-finite 



function of ti, . . . , t k , by Theorem 3. 

We proceed by proving the following invariant of the main loop: the set Ri-± 
generates a D-finite description of 7i(F)(0, . . . , 0, tj, t i+1 , . . . , t k ). This estab- 
lishes the result since it implies that -R^-i contains a D-finite description of 
H(F)(0, . . . , 0, t k ), in this case, a single differential equation. This is precisely 
what the algorithm claims to determine. 

To prove the base case of this invariant, note that Rq contains the generators 
of the intersection (J F W t (t) +X W J n W t (t). We appeal to Corollary 4, to 
conclude that R contains a D-finite description of H(F)(t 1: . . . ,t k ). 



(J F w t (t)+i G )r)W t (t) 



is non-trivial and contains a D-finite description of (F,G). 



Proposition 5 Algorithm 2 terminates and is correct. 



PROOF. First, we remark that for fixed k, 
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The general case is proven with the known result [5] that given a D-finite 
description of a function F(xi,X2, ■ ■ ■ ,x n ), one can compute the D-finite de- 
scription of F(xi, . . . , x n _i, 0), for example, by first eliminating d Xn , removing 
factors of x n in the remaining polynomials, and finally, setting x n = in the 
equations, precisely the process outlined in Algorithm 2. □ 



5 Example: /c-Uniform Tableaux 

Another family of combinatorial objects whose generating function can be 
resolved with our method is a certain class of Young tableaux, namely k- 
uniform Young tableaux. 

For a partition A = (Ai, . . . , A&) h n, a Young tableau of shape A is an array 
T = (T it j) of positive integers T it j defined when 1 < i < k and 1 < j < A«. 
When a Young tableau is strictly increasing on each of its rows and columns 
(Tj j < T i+ ij and T it j < Tj whenever this makes sense) and the n inte- 
gers Tjj are all integers from 1 to n, it is called standard. 

Standard Young tableaux are in direct correspondence with many different 
combinatorial objects. For example, Stanley [34] has studied the link between 
standard tableaux and paths in Young's lattice, the lattice of partitions or- 
dered by inclusion of diagrams. This link was generalized by Gessel [10] to 
tableaux with repeated entries. Gessel remarks that such paths have arisen in 
the work of Sundaram on the combinatorics of representations of symplectic 
groups [36]. 

The weight of a tableau is fi — (/ii, . . . , fik) where fi\ is the number of l's, fi-i is 
the number of 2's, etc., in the tableau entries. Here we consider Young tableaux 
that are column strictly increasing and row weakly increasing, and with weight 
fx = l k 2 k ■ ■ ■ n k : each entry appears k times. We call Young tableaux with these 
properties k-uniform. These correspond to paths in Young's lattice with steps 
of length k. The set of /c-uniform tableaux of size kn is also in bijection with 
symmetric nxn matrices with non-negative integer entries with each row sum 
equal to k. Gessel notes that for fixed k, the generating series of the number 
of /c-uniform tableaux is D-finite [9]. Our method makes this effective. 

Two observations from [21] are essential. First, [xi 1 ■ ■ ■ x% k ]s\ is the number 
of (column strictly increasing, row weakly increasing) tableaux with weight fi. 
Secondly, 

= h[e! + e 2 ] =exp \ Y j p 2 i /2i+ ^ pi/ij , 

A V i i odd / 

which is D-finite. Define yjft to be the number of /c-uniform tableaux of size 
kn, and let Y k be the generating series of these numbers. The previous two 
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1- uniform tableaux 

0o 

01 1 

_<h 

2- uniform tableaux 

00 t 2 {t - 2) 

01 -2(t-l) 2 

J>2 

3- uniform tableaux 

0o (t 11 + 1 10 - 6t 9 - 4i 8 + lit 7 - 15t 6 + 8t 5 - 2t 3 + 12t 2 - 24t - 24) 
<p! -3t(t 10 - 2i 8 + 2t 6 - 6t 5 + 8t 4 + 2t 3 + 8t 2 + 16/ - 8) 

02 9t 3 (-t 2 -2 + t + t 4 ) 

4- uniform tableaux 
4>i (See Appendix A) 

Table 2 



Differential equation 02 F fe " + + 0o^fe = satisfied by fe = 1, . . . ,4. 

observations imply 

n(*) = E# fc = ( ex pfEp?/ 2< + EftA),E fc *-* B ). ( 8 ) 

« \ \i=l i odd / n / 

This problem is well-suited to our methods since again we treat an exponential 
of a polynomial in the p^s, with an explicit closed form in terms of k for this 
polynomial. 

Calculating the equations for k = 1,2,3,4 is fast with either Algorithm 1 or 
Algorithm 2. The resulting differential equations are listed in Table 2. For k = 
1, 2 these results agree with known results [14,34], and are the entries A000085 
and A000985 respectively in Sloane's encyclopedia of integer sequences [32]. 
The first few values of yj^ are summarized in the following table. For k — 3, 4 
these appear to be new. 

Concerning the dual problem, where instead n is fixed and k varies, the se- 
quences (y$^ appear respectively as A019298, A053493, and A053494 for 
n = 3,4,5. Stanley [33, Prop. 4.6.21] reports that the generating functions 
G n (x) = J2k yn^x h are rational with denominator of the form (1 — x) s (l — x 2 ) 1 
where s and t are positive integers. 
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7 M M M 

« Vo ,Vi >2/2 >• 



1 1, 1, 2, 4, 10, 26, 76, 232, 764, 2620, 9496, 35696, 140152, 568504 

1, 1, 3, 11, 56, 348, 2578, 22054, 213798, 2313638, 27627434, 360646314, 

2 5107177312, 77954299144 

1, 1, 4, 23, 214, 2698, 44288, 902962, 22262244, 648446612, 21940389584, 

3 849992734124 

1, 1, 5, 42, 641, 14751, 478711, 20758650, 1158207312, 80758709676, 

4 6877184737416, 701994697409136 



Table 3 

The number, y£\ of fe-uniform tableaux of size kn. 

6 Algorithm for Scalar Product: the General Situation 



So far, we have limited the scope of the algorithms to pairs of D-fmite sym- 
metric functions where only one of the two functions depends on the variables 
ti, . . . , tfc- While this is sufficient in many applications, it is possible to mod- 
ify Algorithm 1 in order to accommodate the t^s in both functions and thus 
make the full power of Theorem 3 effective. While no additional ideas are to 
be used, the description of the algorithm is more technical. 

Algorithm 1 manipulates monomials a and reduces them modulo the ideals 
If and Xq in order to determine equations of the form 

(F, (a - (a red^ T F j) ■ G) = and (f, (a - (a red^ J G )) • G) = 0, (9) 

where on the left, a supposedly does not involve any of the <9 t /s. What makes 
the situation of Algorithm 1 and the left-hand identity in (9) simple is the 
assumption that F does not depend on t, making the action of W t on (F, G) 
act on the right-hand argument only. The difficulty in generalizing lies in that 
now, the action of d ti on F may be non-trivial and must be considered in the 
differentiation rule for scalar products, 

d ti -(F,G) = (d ti .F,G) + (F,d ti .G), (10) 

which itself stems from the differentiation rule for usual products on the level 
of coefficients. 

The idea is therefore to manipulate operators in three sets of dt t 's: one which 
acts on the full scalar product (F,G), and one for each of its components, 
acting directly on the component. To facilitate the description of this situation, 
we denote the former by d ti , the one acting on the left component by dg v and 
the one acting on the right component d n . Using this notation, we wish to 
view Eq. (10) as 

d ti =d ei +d rr (11) 
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We thus modify Algorithm 1 by enlarging the family of monomials over which 
we iterate, and use Eq. (11) to eliminate the <9^'s before we begin Gaussian 
elimination. Here, we iterate over monomials ad^d? of the free commutative 
monoid [p, d p , d^ d r ] with a G [p, d p ] to examine the following generalizations 
of Eq. (9): 

((aftf - (a°3f red&O) • F, d? ■ G) = 

and (d? ■ F, (ad? 

or, with a change of notation, 

(a°dl - (oA9f red£ F ))<9 r 7 ■ (F, G) = 

and fif (ad?. - {ad? red Q G ) ) • (F, G) = . 

Upon making use of Eq. (11) and applying adjunction to the first equation in 
Eq. (12), we get a linear combination of terms of the form <9f • (F, a' ■ G) with 
coefficients in K[t], where 6 N fc , and a' G W Ptt {t). The algorithm proceeds 
as before by performing Gaussian elimination over K(t) to eliminate p, d p , and 
d r . In our implementation, the monomial order ^ is DegRevLex((9 r > de > 
d p > p). The method is summarized in Algorithm 3. 

Algorithm 3 (General Scalar Product) 

Input: F G and G e K[t][[p\], both D-finite in p,t, given by D-finite 

descriptions in W Ptt (t). 

Output: A D-finite description of (F,G) in W t {t). 

(1) Determine a Grobner basis Qg for the left ideal ann^ t ( t ) G with respect 
to any monomial ordering -<, as well as a Grobner basis Qp° for the right 
ideal ann^ p t F° with respect to the same ordering; 

(2) B:={); 

(3) Iterate through each monomial a in p,d p ,de,d r in any order; 

(a) ai := a\d t =dt,d r =i; 

(b) a F := ai - (ai red-< Q F «); 

(c) a r := a\ dr=dti9e=1 ; 

(d) a G := a r — (a r red^ Q G ); 

(e) Introduce (a F \ de=dt - dr )(a\ p=dp=9e=1 ) and (a\ p=dp=dr=1 )a G into B and 
reduce so as to eliminate p, d p , d r ; 

(f) Compute the dimension of the ideal generated by B D W t (t). If this 
dimension is 0, break and output B n W t (t). 

As in Algorithm 1, if m — 1, there is only one variable t, and the condition 
in (3f) is simplified to: 

If B contains a non-zero element P from W t (t), break and return P. 



-(a%Tedg G j)-G) = 0, (12) 
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The same remarks as those made after Algorithm 1 at the end of Section 2 
also apply here. 



7 Termination and Correctness of Algorithms 1 and 3 

7. 1 Sketch of the proof 

The common goal of Algorithms 1 and 3 is to find differential equations satis- 
fied by (F, G), which is equivalent to non-zero elements in W t which annihilate 
(F,G). Although Algorithm 1 is a specialization of Algorithm 3, parts of the 
proof would become artificially more involved if restricted to the simple case. 
We thus treat both algorithms simultaneously. The discussion at the begin- 
ning of Section 2 has illustrated how to manipulate the annihilators of F 
and G to determine a combination P°S + TQ G W t with P G T° F , Q G 1q, 
S G W p (t), T G W Ptt (t), which annihilates (F,G). Not all of the elements in 
ann^ (F, G) are of this form, however, as the following simple example illus- 
trates. If F = pi — p 2 and G = p\ + P2/2, then (F, G) — 1 — 1 = and thus 
1 G annvKt (F, G). However, it can be established that 1 can not be written as 
a combination of the form P°S + TQ for those F and G. Nonetheless, we show 
that the annihilating elements that can be written this way form a non-trivial 
subideal of ann^ (F, G), which we generate with the algorithms. 

Although the problem of finding differential equations appears at first inher- 
ently analytic in nature, we rephrase it algebraically into a question amenable 
to the theory of D-modules. The adjunction properties of the scalar product 
are naturally accommodated by tensor products. Specifically, the proof below 
centers around a certain H^-module S whose elements are tensors, and where, 
for example, 

{i~ x Pi ■ u) <g> v = (u ■ d Pi ) <g> v = u <g> (d Pi ■ v), 

which corresponds to the equivalence ' F,G) = (F,di ■ G). (See also 

Eq. (13-16) below.) On the other hand, the and d n that are involved in the 
description of Algorithm 3 really are the operators d ti <S> 1 and 1 <S> d ti acting 
on S, respectively, where l's denote identity maps. 

The module S can be expressed in terms of the ideal ann Wt (F° <g> G), itself 
contained in ann Wt (F, G). The former ideal is non-trivial and in fact, is suffi- 
cient to describe the scalar product as holonomic, a property whose definition 
is recalled shortly and which implies D-finiteness. In fact, we show that the 
algorithms calculate a Grobner basis for ann Wt (t)(-F <8> G), in other words a 
D-finite description of the scalar product (F,G). 
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The main result is summarized by the following theorem. 

Theorem 6 Suppose F and G are symmetric functions subject to the con- 
ditions of Algorithm 1 (resp. Algorithm 3). Then, Algorithm 1 (resp. Algo- 
rithm 3) determines, in finite time, a Grobner basis for a non-zero D-finite 
ideal contained in annw t (t) (F,G). 

The notion of holonomy to be used in the proof follows [2,7]. Introduce a 
filtration of W t by the fT-vector spaces F d of all operators in W t of total 
degree at most d in t, d t . These spaces are finite-dimensional, of dimension 
(d+2k^ = o(j2k\ ag d tends to i n fi n i t y. a WVmodule M = Y li Wfg i generated 

by a finite family of generators gi is holonomic whenever the .fT-vector spaces 
J2i Fd'Qi have dimension growing like 0(d k ^j . A function of t that is an element 
of a holonomic W^-module is called holonomic. From the definition, it is a basic 
result that a holonomic function is D-finite; the converse is a more difficult 
result to be found in [38, Th. 2.4 and Appendix 6]. Similar definitions apply 
to W^-modules, with a dimension growth of 0{d k+n ^j in place of 0(d k ^j. 

The discussion so far has not relied on the definition of the scalar product. 
Rather, remark that Algorithms 1 and 3 are essentially parameterized by 
the adjunction property of the scalar product of symmetric functions, and 
can easily be redefined and adapted to other adjunctions. It suits our needs 
for the proof to consider adjoints for the usual scalar product of functions, 
{f\d) '■= I f{x)g(x) dx. To avoid confusion, we notationally distinguish (f\g) 
from (F, G) for the two scalar products, as well as * from o for the respective 
adjunction operations. 

Indeed, guided by existing results concerning the preservation of holonomy 
under operations involving the usual scalar product, we link the symmetric 
case to the usual one with a map from one adjunction to the other. This 
reduction also demonstrates how algorithms analogous to Algorithms 1 and 3 
for other scalar products could be shown to terminate with the correct output. 
(See Section 9.1.) 

To make this comparison more intuitive, we could identify (F, G) with the 
integral 

J £(qt-+ F(qi, 2q 2 , ng n )) (p) G(p) dpi • • • dp n , 
where C is the modified Laplace transform 

C(F)(p)= f F(q)e-( piqi+ - +Pnqn) dq, 

which satisfies 

C(q~q i F(qj)(p) = -(a pi o£)(F)(p). 
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Notice, for example: 

(r' P , -F,G) = I C( g „ , ( F( 9l , .... „,„)) (p) G(p) *, ■ ■ ■ 

= - / (d n o £)(F)(p) (a„ ■ G)( P ) dft • • • <i P „ 

= / £(gn F( 9l , . . . , ng n )) (p) (9 W • G)(p) d Pl • • • rfp n = (F, d Pi ■ G) . (13) 

JR™ v 7 

Formally, we must work on the level of abstract modules, however. This avoids 
situations where the integral is not convergent or the Laplace transform is not 
defined as a function. 

Thus, to prove Theorem 6, we show Corollary 10 below which states that 
ann Wt (F° ® G) is a non-zero subideal of axm Wt (F, G) such that the quotient 
W t / ann Wt (F° ® G) is a holonomic module. This is done in several stages. 
First, in Section 7.2, we define S, the algebraic structure in which our cal- 
culations take place, and prove that it is holonomic by reducing the problem 
to the usual scalar product analogue, where similar results are known. This 
analogue is detailed in Section 7.3. Next, in Section 7.4 we express S as a quo- 
tient. Corollary 10 follows from this discussion. Finally, to conclude that the 
algorithm terminates, we relate S to the algorithm in more detail and prove in 
Section 7.5 that all of the generators are determined in finite time. Together, 
these results prove Theorem 6 and thus the correctness and termination of 
Algorithms 1 and 3. 

7.2 The scalar product of symmetric functions 

We now formally define the H^-module S. Begin with U = W Pit ■ F and 
V = W p j ■ G, two holonomic W^t-modules. We shall denote by U the adjoint 
module of U: as K- vector spaces, U = U°, and a right W^^-action is defined 
on U° by u • P = P° ■ u for any u G U° and P G Wp[t], where the last operation 
is taken for the left structure of U. Set S as the tensor product U° ®w p [t] V, 
which makes it a fT[t]-module. This has the desirable effect of encoding the 
scalar product adjunction relations: for all u G U and all v G V, 

(d Pi ■ u) <g) v — (u ■ dp.) (g) v — (u ■ i _1 pi) <g> v = u <g> (i^Pi ■ v), (14) 
(Pi ■ u) <g) v — (u ■ p°) <S> v — (u ■ id Pi ) <S> v — u <g> (id Pi ■ v), (15) 
U • (u <S> v) — (U • u) <S> v — (u • U) <S> v — u <S> (U • v). (16) 

To endow S with a W^-module structure, let d ti act on a pure tensor u <S> v by 
d ti -(u®v) = (d ti -u)®v + u® (d u ■ v), (17) 
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and extend to S by fT-linearity. In other words, = + d Ti after defining 
®ii — d ti <S>1 and d Ti = 1 <8> d ti , where l's are identity maps. 

Armed with this definition and Theorem 7 (formally stated and proven in- 
dependently in Section 7.3), we prove that S is holonomic. Theorem 7 is an 
analogous result for the usual scalar product, corresponding adjunction, and 
corresponding adjoint module M* of a module M. It states that for holonomic 
M and N, M*<S>w p [t] N is a holonomic W^-module under the action of d ti given 
by (17). We shall appeal to this theorem with an appropriate choice for M 
and N. 

To determine the relationship between the two scalar products and make our 
choice for M and N, we compare both adjunction operations. In the sym- 
metric case, adjunction is defined as the anti-automorphism o which maps pi 
to id Pi and d Pi to i~ x Pi-, for all i, and the usual scalar product adjunction is 
defined as the anti-automorphism * which maps d Pi to — d Pi , and leaves the Pi 
variables unchanged. One way to connect both adjunctions is to factor o into 
the composition of three algebra morphisms: 

(1) the automorphism r mapping (pi,di) to (ip i ,i~ 1 d i ). This corresponds to 
the dilation which maps a function F to p F(p 1 , 2p 2 , . . . , np n ); 

(2) the automorphism T mapping (pi,di) to (—di,pi) and named 'Fourier 
transform' in D-module theory (see [2, proof of Th. 3.1.8] or [7, p. 39]). 
Informally speaking, this corresponds to mapping a function F to its 
Laplace transform C(F); 

(3) the anti-automorphism * mapping (pi,di) to (p,, — 9$). 

The important property to note is that each of these three maps preserves 
holonomy since they preserve total degree, hence are filtration-preserving bi- 
jections. A direct calculation on pi and di verifies that o = *ojFor, so that the 
composite o also is a holonomy-preserving linear bijection. Thus, we introduce 
two holonomic modules, M = [T o t)(U) also denoted U Tor , and iV = V, so 
as to appeal to Theorem 7. One concludes that 

S = U* ®w p[ t] V = (U^Y ® Wp[t] V = M* ® Wp[t] N (18) 

is a holonomic W^-module. After we have described the quotient structure of 
S in Section 7.4, this information will be used to prove that annw^-F <8> G) is 
non-trivial and that the quotient module W t / ann^(F° ® G) is holonomic, a 
fact we use to show that the algorithms terminate. 
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7.3 Preservation of holonomy under the usual scalar product 

In the previous section, we reduced the proof of the holonomy of S — U°®w p [t\ 
V to an analogous result in terms of the usual scalar product, to be proven in 
this section: the module T = M* ®w p [t] N is holonomic when M and iV are. 

The following notion will be used in the proof: the integral of a W^-module P, 
denoted J P = J P dpi ■ ■ -dp n , is defined as P j (j2id Pi ' P)- It is the im- 
age of composed maps: the Fourier transform J 7 , the inverse image 7T* under 
the projection ir from K n+m to K n defined by n(p,t) = t, and the inverse 
Fourier transform. Specifically we have, / P = J r ~ 1 7r*J-'(P). These maps pre- 
serve holonomy (see [2, Th. 3.3.4] or [7, Th. 18.2.2 and Sec. 20.3]), so that the 
integral of a holonomic W P:t -modu\e is a holonomic W t -modu\e. (See also [2, 
Th. 3.1.8].) 

The module T fits naturally in between an existing holonomy-preserving sur- 
jection from the WVmodule / M ®K]p,t] N to the space (M\N). Factoring this 
map to pass through T = M* ®w p [t] N yields: 

J M ® KM N^UM* ® Wp[t] N (M\N) , (19) 

where ip surjectively maps m ® n to (m\n), and is a natural W^-linear 
surjection that we are about to define in the course of the next theorem. After 
proving that the first module in (19) is holonomic, the surjectivity of (f> implies 
the holonomy of T. 

Theorem 7 Suppose that M and N are two holonomic W Ptt -modules, and 
define T as M* ®w p [t\ N. Then, T is a holonomic W t -module under the action 
of d ti given by 

d ti • (m <S> n) — (d ti • m) <S> n + m <S> (d ti ■ n). 

PROOF. First, we focus our attention on the module / M ®K\p,t] N m (19). 
Consider the W^, it -module P := M ®K\p,t] N, with action of d Pi defined by 
d Pi ■ (m <S> n) = (d Pi ■ m) ®n + m® (d Pi ■ n), and action of d ti defined similarly. 
We can also write this as the inverse image i* (M ®x N), where i is the map 
from K m+n to f^(™+ m )+(™+ m ) which sends (p,t) to (p,t,p,t). The advantage 
of the second presentation is that the holonomy of P is obtained from the 
holonomic closure under inverse image under embeddings (see [2, Th. 3.2.3] or 
[7, Sec. 15.3 and Ex. 15.4.5]) and the holonomic closure under tensor product 
over K [7, Cor. 13.4.2]. Therefore, f P is also holonomic. 

Next, we define a PVt-linear surjection to T. Define a map from M x N to T 
which sends (m, n) to m®n. This map is K\p, t]-balanced, K[p, t]-bilinear, and 



31 



surjective. By the universality of the tensor product, this induces a surjective 
map from P = M <8>K\p,t] N to T. Observe that each derivation d Pi maps P 
into the kernel of 0, as the following calculation indicates: 

(f)(d Pi ■ (m <g> n)) = (j)({d Pi ■ m) ® n + m ® (d Pi ■ nfj 

= (d Pi ■ m) ®n + m® (d Pi • n) — m<g> (—d Pi ■ n) + m® (d Pi ■ n) = 0. 

In other words, J2i 9 Pi -P C ker 0, and thus also induces a well-defined surjec- 
tive map from / P to T. Any good filtration of / P will induce a good filtration 
for T (see [2, Prop. 1.11] or [7, Lemma 7.5.1]). Thus, T is finitely generated 
with dimension bounded by that of / P. Therefore, T is holonomic. □ 



7.4 The quotient structure of S 

Subsequent developments to express S as a quotient involve modules over W Ptt 
and ideals of W p j, rather than W p j(t). We therefore introduce the annihi- 
lators If = ann^ ( F and Iq = ann^ t G, to be used in place of Xp = 
&nn Wp ( t ) F and Xq = ann^^G, respectively. Note that Ip — Xp n W Pt t 
and Xp = K(t) <S>kh] If, and similarly for G. Finally, although adjunction has 
not been defined for d t , we use the notation W^ t to denote W Ptt endowed with 
both a structure of Wt-module on the left and a structure of W 7 ^] -module 
on the right. 

Proposition 8 The module S = (W Pjt ■ F)° ®w p [t] (Wp,t ' G) is isomorphic to 
(W p t <S>w p [t] W p t )/(I F ®w p [t] W p t + W p t 

®w p [t] Ig)- 



PROOF. The W t -module S = U° <S> Wp [t] V is also a W* t ®w p[t] VF^-module. 
As such, it is generated by F° <g> G. Consider the two exact sequences of 
respectively right and left W p [t] -modules 

-> I F A W« t A U* -> 0, 
- i G A W Ptt ^ V - 0, 

where ot(P) — F° ■ P, (3{Q) = Q ■ G, and p and a are inclusions. (Here, F 
and F° denote the same element of the set U, but we write F° when viewed as 
an element of the right module U°, F when viewed as in the left module U.) 
We combine them to make a third exact sequence: 

ker(« ® 0) - W; tt ® Wp[t] W ptt ^ S -0, 

P®Q .— > {F« ■ P) ® {Q ■ G) 
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where, by [3, 11.59, Proposition 6], 

ker(a <g> (5) = im(p <g> l Wpt ) + im(l w * t <g> a) = I F ® Wp[t] W pt + W* t ® Wp[t] I G 
as K[t] -modules. We conclude that, as W t -modules, 
S-(W« t ® Wp[t] W pit )/ker(a®P) 

* (w° t ®w v[ t] w p ,)/(i F ® Wp[t] w p , + w; tt ® Wp[t] i G ). 

□ 

To be more explicit, note that this isomorphism maps the class of 1 <8> 1 in the 
quotient to F° <g> G G S. Remark also that, as W^-modules, 

ker(a ®/5) = {?®Qe W* t <g> W p , t ■ (a <g> (3){P ® Q) = o} 

= {p®Qe W; <t <g> W p>t : (F° ■ P) <g> (Q • G) = 0} 
= {P ® Q G Wp° >t ® Wp, t : (P (8) Q) ■ (F° ® G) = o} 

= anii W o t8wi>[(]Wj>it (i^(8)G), 

so that we also have 

aim^ t ® Wp[t] ^ (21) 

Proposition 9 The W t -module S' — W t ■ (F° ® G) is a submodule of S, 
isomorphic to 

K I ((If ®w p[ t] W p:t + W° t ® Wp[t] I G ) n W[) , 

where W[ ~ £/ie smallest K-subalgebra of W p j <S>w p [t] W Ptt generated by 

K[t], 1 ® d tl + d tl ® 1 , ■ ■ ■ , 1 ® <9 tfc + 9 tfc ® 1 . In the simplified situation when 
If = d t W p>t + W t Jp for J F = &nn Wp F, S' is isomorphic to 

w t I ((w t jp + i G )nw t ). 

We first prove this proposition, then in the next section we discuss how to 
connect the description of S' above directly to the algorithm and how to 
apply it to show that the algorithms terminate. 

PROOF. The annihilator of F° ® G in W( ■ (F° ® G) 

ann^(P° ® G) = ann^^^F ® G) H W[. 
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In view of the action of Wt on S' through the isomorphism between Wt and W{, 
we thus have that S' is isomorphic to W t / ann W( (F° <S> G) , itself isomorphic to 

Wl/ajm w ,(F* ® G) = W[j ' (&nn w « t ® Wp[t]W jF« ® G) n w{). 

Owing to (21), this proves the general quotient expression for S' in the propo- 
sition statement. 

Now, to prove the formula in the simpler case, observe that when I F = d t W p>t + 
W t J F , 

If ®w m w Ptt = d t w; tt ® Wp[t] w Ptt + w t r F ® Wp[t] w Ptt 

= d t w t ® K [t] w Pit + w t ® K [t] w t r F 

while W£ jt ® Wp [ t ]I G = W t <8) K [ t ]I G , whence the relation ker(a(g)/3) = d t W t <S>K[t] 
W Ptt + W t ® K[ t] (W t J° F + I G ). Since W£ t ® Wp [t] W p>t = W t ® K[t] W Ptt , we obtain 

S~ W p , t /(W t J F + I G ), 

as (W t ® K[t] W Ptt )/kev(a®P) ~ (K[t\ ® K[t] W Ptt )/(K[t] ® m (W t J F + I G )) ~ 
W Pit /(W t J F + I G ). Following these isomorphisms, W[ can be identified as the 
copy of W t included in W Pit in the last quotient above. Therefore, the sub- 
module 5" of S is isomorphic to the quotient announced in the proposition 
statement. □ 

Corollary 10 The ideal ann Wt (F° (g) G) is: 

(1) isomorphic to (I F <S>w p [t] W p t + W° t <S>w p [t] Ig) H W' t as a W t -module; 

(2) a non-trivial ideal contained in ann^ ( (F,G) and such that the quotient 
W t / anmy^-F <8> G) ~ S' is holonomic. 

PROOF. From (21), 

ann^(F° ® G) = (ann^ t ^ [t] ^(F° <g> G)) n W[ 

= (l F ®w p[t] W p , + W; :t ® Wp[t] l a ) n Wl, (22) 

and we have shown (1) in the corollary statement. The WVmodule S' ~ 
W t / ann^F <g> G) is a holonomic WV module, as it is a sub module of the 
holonomic W^-module S. Now since W t is not holonomic, ann^ (F° <g) G) must 
be non-trivial by a simple dimension argument. Finally, we recall that this non- 
trivial ideal is contained in ann^ (F, G), since there is a surjection from S' to 
W t / annyKi (F, G) given by ip : (u®v) h- > (u, v). This proves (2) in the corollary 
statement. □ 
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7. 5 Termination 



We now link the modules S and S' to the algorithms and prove their termina- 
tion. The termination of Algorithm 3 is more technical to prove than that of 
Algorithm 1 since can act separately on F and G. Thus, for ease of presen- 
tation, we consider Algorithms 1 and 3 in turn, to show that they eventually 
generate a Grobner basis for ann Wt ^(F° <g> G). 

7.5.1 Termination of Algorithm 1 

The basic idea of Algorithm 1 is to compute nitrations of T F and T G indepen- 
dently and incrementally and to recombine them at each step. The algorithm 
terminates when condition (3e) in the algorithm description is satisfied. We 
show that the algorithm will satisfy this condition by eventually producing a 
Grobner basis for ann^^F* <8> G). This subideal describes F°®G and (F, G) 
as D-finite. 



PROOF. (Theorem 6, Algorithm 1) Algorithm 1 places a constraint on F 
that allows us to take advantage of the simpler ^-structure of U — W Pjt ■ F: 
since each d ti ■ F is 0, we have U = K[t] ®k (W p ■ F) and I F = d t W p>t + 
W t J F . Taking the intersection with W[ is then far more transparent: from the 
previous section, we obtain the following simplification of Eq. (22): 



Considering the monoid of monomials generated by p, d p , d t , ordered by the 
monomial order ^ specified by the algorithm, we denote by Vp the filtration 



Assume that Algorithm 1 fails to terminate on some input F and G. For any f3, 
Algorithm 1 thus eventually reaches a value for the main loop index a such 
that all the monomials that have been considered in the algorithm span a 
vector space containing Vp. After Step (3d) in the main loop for this value a 
of the loop index, B generates a vector space containing 



By our choice of elimination term order, B fl W t (t) consists of generators of a 
vector space which contains the intersection Lp fl W t (t). 



Next, for each 7, ( J F W t (t) + TgJ fl V 7 is a subspace of Lp for some (5. In- 
deed, since V 7 is finite-dimensional, so is the intersection under consider- 
ation. Let us introduce a basis bi,...,b d of it. Each 6j can be written in 



ann m (F° ®G) = (J F W t + I G ) n W t . 



(23) 



Lp := (r F w t (t) n v^) + (la n Vp). 
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the form /j + g,i for f\ & T F = J F W t (t) and ^ G X G , so that, provided 
P = max{maxj deg /j, maxj deggj}, the intersection 

(j£w*(t) +z G ) n v 7 = ©tf(t)(/i + 

i=i 

is a subspace of 

d d 

Y,K(t)fi + Y,K(t) 9l c (w t (t)r F n ty) + (x G n ty) = v 

i=l i=l 

Since ann^^F ® G) is finitely generated by noetherianity of Wt(t), we can 
choose a finite set of generators for it, and set 7 to their maximal leading 
monomial. Consequently, the chosen generators are in 

ann m(t) (F <> <g> G) n V 7 = (w t (t) J£ + X G ) n W t (t) n V r 

By the reasoning above, the latter is a subspace of Lg for some (3, and 
when the loop index reaches a sufficiently high a, annw t ( t ) (F* <8> G) is a 
subideal of the ideal generated in W t (t) by B D Wt(i). Since, by Corollary 10, 
W^/ ann^ (F° <8> G) is a holonomic module, ann^^F ®^) is of dimension 0, 
and condition (3e) is satisfied. The algorithm terminates, a contradiction to 
our assumption. □ 

A limitation of the algorithm is that we cannot predict in advance how many 
monomials must be tested, and hence cannot estimate the running time. 

7.5.2 Termination of Algorithm 3 

The termination of Algorithm 3 can be proved similarly, but we must use 
greater care when treating the d ti . 

PROOF. (Theorem 6, Algorithm 3) Since there is no adjoint action for dt t , we 
consider occurrences of in the left argument of the scalar product differently 
from those on the right side. This is modelled in S by tensoring over Hp[i], 
where d t is absent and thus, d ti ®l differs from l<S>d ti . Both still obey the same 
commutation law with tj as d ti . Denote the former by and the latter by d n . 

Having distinguished these two cases, we rewrite several of the important 
elements from the previous proof using this new notation. For example, 

W£ t ®w p [t] W pt = K(p,t,d p ,d e ,d r ; [d Pi , Pj ] = [d t .,tj] = [d r .,tj] = 5 id , 

[PhPj] = \Pi,tj] = [U,tj] = [d e .,pj] = [d r .,pj] = [d Pi ,tj] = 0), 
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d ( 



-d rk . We can also 



and its subalgebra W[ is generated by K[t], de 1 + d, t 

rewrite 1% ® Wp[t] W p>t + W* t ® Wp[t] I G in the form P F dt=d K[d r ] +K[d e ]I G 
Algorithm 3 actually computes with coefficients that are rational functions 
in t, and so with elements of X F K[d r ] + K[dg\X G 



dt=d r 



dt=d e 



dt=d r 



In order to endow W£ t ®w p [t] W t with a filtration, let us extend the ordering ^ 
to monomials in p, d p ,de, d r by considering any ordering which, after setting 
de = dt, d r = 1 or d r = dt, de = 1, respectively, induces the ordering -<. We 
denote the extended ordering by -< as well. Then, we let Up denote the filtration 
7 ^ K{t)f3 for j3, 7 ranging over the monomials in the variables p, d p , d r , de- 
Turning our attention to W((t), let be the image of the Vp of the previous 
section, under the same transformation which takes W t (t) to Wl(t), that is, 



K(t)p a d b p (d e + d r y 



p a d p d^p 



For each f3, there is (5' such that V'p C U l 



Assume that Algorithm 3 fails to terminate on some input F and G. Since the 
main loop enumerates all monomials in p, d p , de, d r in some order, for any f3 
there exists a value of the index loop a such that when the loop reaches it, 
all monomials that have been enumerated span a vector space containing Up. 
After the algorithm has introduced (variants of) ap and a G at Step (3e) for 
this value of a, let us call V a the vector space generated by the set B. Setting 
de = d t — d r maps V a to a vector space which contains 







dt=d t 



K[d r })nUp+(K[d e ]X G 



dt=d r 



nu 



■(3- 



We use this fact to conclude termination. 



At this point we show that for each 7, the vector space X D where 

x = xp ® Wp{t) w Ptt (t) + w Ptt (ty ® Wp{t) x G 

is a subspace of Hp for some f3. Indeed, choose 7' such that Vi f C Uy, so 
that X fl Viy C X DUji. The latter intersection is finite-dimensional, since 
Uy is so. Suppose it has for basis bi, . . . ,bd, with each bi of the form b; L = 
fiU + Ugi, where fiEX F g { E X G , n G K[d r ], and l L G K[d t ], and 

d t =Ot dt=d r 

set (5 = maxjmaxj deg f^i, maxj deg/j^}, where here deg extracts the leading 
monomial. Then, 



xnv'pc® K(t)(f iri + l i9i ) c K(t)Un + £ K(t)l i9i c Hp 



i=i 



i=i 



i=i 



By noetherianity, we can choose a finite set of generators for ann W / i ( t )(F°®G), 
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and set 7 to their maximal leading monomial. The generators are thus elements 
of ann m(t) (F°(g)G')nV 7 , which is isomorphic to &mi w ^ t) (F ®G) H V' y . By (22) 
the latter is also X n V 7 , and, as explained above, there is (5 such that this is 
a subspace of Hp. 

By our earlier loop invariant, the same generators, after setting dp = d t —d r , are 
contained in the space spanned by B when the loop index reaches a sufficiently 
high a'. Thus, it suffices to run the algorithm until this a and generators of 
&nnw t (F°®G) will be contained in B. At this point the termination conditions 
are satisfied, and the algorithm terminates. □ 



8 Asymptotic Estimates 

We now illustrate how the differential equations computed by our algorithms 
may be exploited in order to derive asymptotic estimates of combinatorial 
quantities. 

8.1 Outline of the method 

A very general principle in asymptotic analysis is that the asymptotic be- 
haviour of a sequence is governed by the local behavior of its generating series 
at its singularity of smallest modulus, see for instance [28, Section 10]. Our 
approach is thus based on applying the classical analysis of linear differen- 
tial equations as presented in textbooks such as [17,43] in order to derive 
asymptotic estimates for the coefficients. Moreover, large parts of this analy- 
sis can be automated thanks to the algorithms described in [23,39,42], many of 
which have been implemented in computer algebra systems 6 . An alternative 
approach based on Birkhoff's work can be found in [44]. 

In the special case of solutions of linear differential equations, the possible 
location of singularities is restricted to the roots of the coefficient of the high- 
est derivative. Then, the analysis depends on the nature of the singularity. 
The classical theory distinguishes two kinds of singular points: regular sin- 
gular points, where the solutions have an algebraic-logarithmic behavior; and 
irregular singular points where the solutions have an essential singularity of 
the type exponential of a rational power. Accordingly, the asymptotic behav- 
ior of the coefficients is deduced either by singularity analysis [8,18], or by the 
saddle-point method [16,45]; both approaches are implemented in the algolib 
library. 

6 In Maple, this functionality is provided by DEtools [f ormal_sol] . 
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This asymptotic analysis of D-finite generating series extends to the divergent 
case. Indeed, the coefficients u n of a divergent D-finite series grow at most like a 
power of n\ with a rational exponent p/q which can be computed (see example 
below). Then one constructs an auxiliary differential equation satisfied by the 
convergent generating series of u n /(n(n — q)(n — 2q) ■ ■ -r) p (where r denotes 
the remainder of the division of n by q), to which the previous method applies. 
This construction is achieved thanks to the closure properties of D-finite series, 
by multiplying u n with the solution of the recurrence (n + q) p v n+q = v n , which, 
up to a constant, grows like n\ p / q n p ( q ~ 1 ^ 2q . This operation is implemented in 
the gf un package. 

8.2 k-uniform Young tableaux 

We now illustrate this method in the special case of the fc-uniform Young 
tableaux of Section 5. We treat in detail the case k = 3; other cases are 
similar. To the best of our knowledge, these asymptotic estimates are new. 

We start from the differential equation for k = 3 to be found in Table 2. This 
is a second-order differential equation and its leading coefficient vanishes at 
the origin. This indicates a possible singularity of Ys(t) at the origin, which 
would be reflected by the divergence of this series. Indeed, from this differential 
equation, a linear recurrence is readily computed for the coefficients u n := yjf': 

u n + u n+ i - (3n + 12)w„ +2 - 4w n+3 + (6n + 35)w n+4 - 15w n+5 
+ (9n 2 + 93n + 242)w„+ 6 + (18n + 126)w n+7 - (9n 2 + 159n + 698)w n+8 
+ {9n 2 + U7n + 606)w„ +9 - (18n 2 + 366n + 1884)w n+ i 

- (48n + 552)w n+11 + (24n + 288)w n+ i 2 = 0. 

8.2.1 Divergence 

From this recurrence it is easy to compute a couple hundred coefficients and 
observe their rapid growth. Simple experiments indicate that the growth of 
these coefficients is of order \fn}.. That this growth is the exact exponent of n\ 
in the behavior follows upon considering the degrees of the coefficients in the 
recurrence: the terms of order 12 and 11 have coefficients of degree 1, while 
the term of order 10 has a coefficient of degree 2 (the maximal degree). Thus, 
up to first order, the behavior is dictated by 

24nw„ + i2 = 18n 2 M„ + i , 

which leads to a growth of order (^) n / 2 n\ 1 / 2 . In order to derive a more precise 
estimate, we compute a linear differential equation satisfied by the convergent 
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generating function of y$v n where v n satisfies v n+2 = v n / (n + 2). This differ- 
ential equation is obtained by first computing a linear recurrence for y$v n , 
which exists thanks to the closure properties of linear recurrent sequences. This 
closure operation produces a linear recurrence of order 24 with coefficients of 
degree 29. From there we obtain a linear differential equation of order 29 with 
coefficients of degree 37, which we now analyze. 

8.2.2 Singular behavior 

The leading coefficient of the previous equation is t 27 (3t 2 — 4), up to a constant 
factor. This reveals a dominant singularity at p = 2/^/3, thus confirming 
the growth order (3/4) n//2 expected from the previous stage 7 . The next step 
consists in analyzing the behavior of our convergent generating series in the 
neighborhood of p. A local analysis of the differential equation reveals that all 
solutions of this equation of order 29 behave like 

, x ^ ex P(ll) /, 145 8591 2 ox\ , 

giu) + i/u~ ( - 144" - 41472" + ° (" )) ' 1 - z/ " = " - °' 

where g is an analytic function at 0, and A is a constant depending on the 
solution. 

8.2.3 Asymptotic estimate 

This behavior is typical of an irregular singular point and can thus be dealt 
with using the saddle-point method. Putting everything together, we finally 
obtain 

for some constant C3, and where the O-term hides the beginning of an ex- 
pansion in descending powers of n that could be computed with the same 
method. 

The constant C3 can then be approximated numerically by using Romberg's 
acceleration method, adapted to powers of n -1 / 2 , and we get: 

C 3 w 0.377200. 

8.2.4 Other values of k 

The computation of the asymptotic behavior of yj® for other values of k is 
completely similar, provided one has computed the differential equation. We 

7 We could also have incorporated this factor in the recurrence for v n . 
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1 c x d« 0.347829 

vn! n 1 / 4 



2 C 2 eXP )_ 2n C 2 « 0.282094 



0.377200 



o n /—,( ex P ^ n 
3 <W^-j C 3 



( 2V'oxp^7 r|%()iS . 5| - (i - 



Table 4 

Asymptotic number of /c-uniform Young tableaux 



summarize our results in Table 4. This serves to illustrate a typical use of 
our techniques in experimental mathematics to obtain conjectures such as the 
following. 

Conjecture 11 The number y$ of k-uniform Young tableaux of size n be- 
haves asymptotically according to 



This conjecture is proved for k — 1 and k = 2: the constant is obtained from 
a closed form solution of the differential equation. For k = 3 and k — 4, only 
the value of the constant is conjectural. The proof of the general case of the 
conjecture requires techniques such as those of [11,25], which fall outside of 
the scope of this article. 



8.3 Conclusion 



The main advantages of our method are its general applicability, its ability to 
produce full asymptotic expansions up to one constant factor, the availability 
of computer algebra programs that automate many of its steps. The price to 
pay for this generality is that the method can only produce numerical estimates 
for the constant factor. In some special cases, specific approaches often exist 
that provide this constant term. 
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9 Conclusions and Directions for Future Work 

9. 1 Applying the method to other scalar products 

Let us note that the method of this article can be applied in the case of 
other scalar products, provided that the corresponding adjunction o (no longer 
denoting the symmetric adjunction) is a linear involution that preserves the 
total degree (in p, d p ) of the differential operators. In effect, one should simply 
set M = (U°Y and N = V to obtain a suitable analogue to (18) and prove the 
holonomy, thus D-fmiteness, of the scalar product: M is holonomic if and only 
if U is. Since the statement and proof of Algorithm 1 and 3 do not make use 
of any other special property of o than being a degree-preserving involution, 
correctness of the algorithms can then be established along the same lines as 
for the case of the scalar product of symmetric functions. 

We use this idea in the next two sections by introducing various scalar products 
given by an adjunction relation involving a formal parameter. 

9.2 Calculating the Kronecker product of symmetric functions 

Another symmetric function operation, closely related to the scalar product, 
is the Kronecker product, also known as the tensor product. One can define 
it on the power basis as p\*p ll = (p\,Pn)p\- Gessel showed in [9] that given 
two D-finite symmetric series F and G, the Kronecker product F * G is also a 
D-finite symmetric series. Algorithm 1 can be used to make this fact effective 
via the following observation: 



More precisely, we rewrite a Kronecker product as a scalar product by multi- 
plying each pi in F by ij. In the system which results we make the substitution 

U = pi and d u = d Pi . 

We formalize this in the following algorithm, which merely calls Algorithm 1 
on modified input systems. 

Algorithm 4 (Kronecker Product) 

Input: Symmetric functions F G if [[p]] and G G if [[p]], both D-finite in p, 
each given by a D-finite description in W p . 
Output: A D-finite description of F *G in W t . 

(1) Call Q the system defining G and set Q' = {tid tl — pid pi , . . . , t n d tn — 
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Pnd Pn }; 

(a) For each element in Q , replace pi with Upi, d Pi with t~ 1 d Pi and add 
toQ'; 

(b) For each element in Q, replace Pi with t(Pi, d Pi with p~ x d ti , clear 
denominators, and add to Q 1 ; 

(2) Follow the steps of Algorithm 1 on the input system for F and the modified 
system Q' for G; 

(3) In the output of Algorithm 1 make the substitution ti = Pi and d ti = d Pi 
and return this value. 

Many interesting problems which use this operation require an infinite number 
of p n , and are thus at first glance seemingly unsuitable for direct application 
of our algorithms. However, applying our algorithms for several truncations of 
a combinatorial problem can serve as a means to generate information upon 
which reasonable conjectures can be formulated. For example, Eq. (25) below 
was initially conjectured after a clear pattern emerged from a sequence of ap- 
peals to Algorithm 4. For each of these, we render the problem applicable by 
setting most p n 's to 0. In some cases, notably symmetric series arising from 
plethysms, there is sufficient symmetry and structure which can be exploited 
to verify these guesses by applying one of Algorithm 4 to well chosen subprob- 
lems. That is, in certain cases, such as the example that follows, the Kronecker 
product of two functions each with an infinite number of p n variables can be 
reduced to a finite number of symbolic calculations. 

For example, if two symmetric series F and G can be expressed respectively 
in the form 

F(pi,p 2 ,- ■ ■) = Yl fn(Pn) and G(p u p 2 , . . .) = JJ 9n(Pn), 

n>l n>l 

for functions /„, g n , then one can easily deduce that 

F*G= Uf n (Pn)*9n(Pn). (24) 
n>l 

Remark that series which arise as plethyms of the form h[u] or e[u], where u 
can be written as a sum J2n u n(p n ), for some functions u n , are precisely of this 
form. For example, we can use this fact to compute the Kronecker product of 
the sum of all Schur functions 

F( Pl ,p 2 , ■ ■ •) = E *a = %i + V2p? - l/2p 2 ] = exp (E | + ^ 



2i 2i - 1 / ' 



and itself. Due to the patterns present, we can reduce the calculation of the 
entire product to two symbolic calculations. More precisely, in order to deter- 
mine a system of differential equations satisfied by G = F * F we consider 
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only the even and odd cases, and set 

f 2n = exp(p 2 n /4n) and f 2n -i = exp((p 2 2n _ l /2 + p 2n -i)/ (2n - 1)). 

All of the functions #2n = jin * fzn are obtained from a single computation 
by our Algorithm 4, adapted to handle a formal parameter. This modification 
is of the same nature of that described in Section 9.1. Here we introduce the 
scalar product given by the adjunction formula p° = nd for a formal parameter 
n from the field K. Thus computing exp(p 2 /4n) *exp(p 2 /4n) with this variant 
algorithm results in a first-order operator in p and d, which, once interpreted 
back in terms of p n becomes: 

(1 - P n )—K 1- Pn9n{Pn) = 0, for even n. 

Op n 

A second calculation for g 2n -\ = fin-\ * iWi-i results in: 

n{\ +p n )(l -p n f^PA - (l + ( n + l) Pn - np 2 n ) g n (p n ) = 0, for odd n. 

OPn 

These linear equations are satisfied respectively by the functions 

-1/2 



92n ={}- P2n) aIld 

V2n-\ \ (, 2 \-V2 

-)) (i-ftn-lj 



^2n-l = CXp 



(2n-l)(l-p 2f 
Applying Eq. (24) above, we get the following result. 

Proposition 12 TTie Kronecker product of the sum of the Schur functions 
with itself is 



-1/2 



(25) 



9.3 A q-analogue 



A g-calculus parameter can be incorporated in symmetric functions in sev- 
eral ways. 

Apart from the scalar product defined by (1), several other ones are of interest 
in relation to symmetric functions, notably the following two, which lead to 
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the definitions of Hall and Macdonald polynomials respectively: 

(Pl»P\) =2A<Wn( 1- * Ai ) alld (Pm^Pa) = Z X S^x II i \ i i 

i=i i=i 1 9 

where £(A) is the length fc of a partition A = (Ai, . . . , A&). The same approach 
as in this article works in this setting and our Maple code has been adapted 
very easily 8 . 

As a related problem, the ring homomorphism 9 q : A — > if [?][[*]] defined as 
e q (f( Xl , x 2 , ...))= /((I - q)t, (1 - q)qt, (1 - q)q\ . . .) 

is useful for studying partitions and for counting permutations [34]. This is 
one possibility for a g-analogue to the map 9 from Theorem 1 (named expo- 
nential specialization in [34]), since \\m. q ^ 1 9 q (F) = 9{F)(x). An algorithm to 
compute 9 q , possibly mapping differential equation to D q equation should be 
of interest. 

9.4 Other conditions for D- finite closure 

Remark that Theorem 3 requires that g be a function of only a finite number 
of p n . The necessity of this condition is evident in the following example. Find 
a sequence c n such that J2 c n t n is not D-finite. However, according to the given 
definition of D-finite symmetric series, J2 n c nPn is D-finite, as is J2 n Pnt n /n. The 
series (Y,n c nPn,Y, n Pnt n /n) = J2n c nt n is not D-finite by construction. 

On the other hand, the condition is not essential. We have that (H(l), H(t)) = 
which is D-finite despite H being a function of all p n . Perhaps a closer 
investigation on the level of modules could reveal a refined condition. 
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A 4-Uniform Young Tableaux 



The differential equation satisfied by Y A (t) is 

64t 4 (t - 2) 2 (t + l) A a(t)Y^(t) - lQt 2 (t -2)(t+ l) 2 /3(t)F 4 (2) (t) 

+ 4rf(t)Yi(t)-S(t)Y A (t) = 

where a(t), f3(t), j(t), S(t) are irreducible polynomials given by 

a(t) = t u - t 13 - 5t 12 - 7t n + 6t 10 + 35t 9 + 39t 7 - 50t 6 - 162t 5 - 92t 4 

+ 228t 3 + 424t 2 + 248t + 48, 
/3(t) = t 29 - 3t 28 - 16t 27 + 24t 26 + 147t 25 + 14t 24 - 770t 23 - 666t 22 + 1416t 21 

+ 3567t 20 - 916t 19 - 16598t 18 + 17766t 17 + 40678t 16 - 102556t 15 

- 53272t 14 + 390656t 13 + 364080t 12 - 707936t n - 1406336t 10 - 552544t 9 
+ 1397664t 8 + 2020864t 7 + 176256t 6 - 916864t 5 + 304896t 4 + 1283328t 3 
+ 877056t 2 + 253440t + 27648, 

7 ( t ) = t 28 - t 27 - Ut 26 - 20t 25 + lilt 24 + 278t 23 - 196t 22 - 1216t 21 

- 1384t 20 + 2765t 19 + 3170t 18 - 3400t 17 + 12140t 16 + 15588t 15 

- 70280t 14 - 108946t 13 + 121796t 12 + 349056t n + 116992t 10 - 481704t 9 

- 706320t 8 + 3040t 7 + 581184t 6 + 158688t 5 - 297408t 4 - 173952t 3 
+ 22272t 2 + 35712t + 6912, 

6(t) = 2t 21 - 3t 20 - 17t 19 - 2t 18 + 74t 17 + 105t 16 - 108t 15 - 172t 14 - 252t 13 
+ 432t 12 - 667t n + 1500t 10 + 7336t 9 - 3772t 8 - 23056t 7 - 20584t 6 
+ 15504t 5 + 38160t 4 + 17904t 3 - 4512t 2 - 5568t - 1152. 



B Sample Maple Session for 3-Regular Graph Computation 



The following Maple session indicates the user-level routines required to pro- 
gram Algorithm 2. It requires the library algolib, which is available at 
http : //algo . inria.fr/packages/. 

# Load the packages . 

with(Ore_algebra) : with(Mgfun) : with (Groebner) : 

# Determine the DE satisfied by the generating function 

# for 3-regular graphs. 

k:=3: Fp:= exp(l/2*pl~2-l/4*p2~2-l/2*p2+p3~2/6) : 
Gp:=exp(l/6*t3*pl~3+l/2*t2*pl~2+tl*pl+l/2*t3*p2*pl 
+l/2*t2*p2+l/3*t3*p3) : 

# Define the variables. 

vars:= seq(p||i, i=l..k): dvars:= seq(d||i, i=l..k): 
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tvars:= seq(t||i, i=l..k): dtvars:= seq(dt||i, i=l..k): 

# Define the algebra. 

A:= dif f _algebra(seq( [dvars [i] , vars[i]], i=l..k), 
seq( [dtvars [i] , tvars [i] ] , i=l . . k) , polynom={vars}) : 
At:= dif f _algebra(seq( [dtvars [i] , tvars [i]], i=l..k)): 

# Define the monomial orders. 

T[g] :=termorder(A, lexdeg( [dvars , vars] , [dtvars] )) : 
T[f ] :=termorder(A,tdeg(vars, dvars, dtvars)): 

# Define the systems. 

sys [g] :=df inite_expr_to_sys(Gp, F(seq(p I I i : : dif f , i=l..k), 

seq(t| |i: :diff , i=l. .k))) : 
newsys [g] :=subs( 

[seq (dif f (F( vars , tvars) , vars [i] )=dvars [i] , i=l . .k) , 
seq(dif f (F(vars , tvars), tvars [i] )=dtvars [i] , i=l..k), 
F (vars, tvars) = 1] , sys [g] ) : 

# Find the Groebner basis for G. 
GB [g] : =gbasis (newsys [g] , T [g] ) ; 

# Do the same for F. 

sys [f ] :=df inite_expr_to_sys(Fp, F(seq(p I I i : : dif f , i=l..k))): 
newsys [f ] : =subs( [seq(dif f (F(vars) , vars [i] )=dvars [i] , i=l . .k) , 
F(vars)=l] ,sys [f ] ) ; 
GB [f ] : =gbas is (newsys [f ] , T[f ] ) ; 

# Define the adjoint and reduction procedures. 
star:= x->subs( 

[seq(d| |i=l/i*p| |i, i=l. .k) ,seq(p| |i=d| |i*i, i=l..k)],x): 
rdc[f] := x->star (star (x) -map (normalf , star(x) , GB[f], T[f])); 
rdc[g] := x->normalf(x, GB [g] , T[g]); 

# Reduce the Groebner basis of F. 

for pol in GB[f] do m[pol] :=rdc[g] (pol) end do: 

# Small optimization: we will always try to reduce with respect 

# to a linear term when possible. 

lpol: = [seq(m[i] , i=subsop(l=NULL,GB [f ] ) ) ,m[GB [f ] [1]]] : 

for indelim from k-1 by -1 to 1 do 
# eliminate dt . indelim 
for j from 2 to nops(lpol) do 

newpol [j] :=skew_elim(lpol [j] , lpol [1] , dt I I indelim, At) 
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end do; 

# set t . indelim = 

lpol : =map (pr impart , subs (t I I indelim=0 , 

[seq(newpol [j] , j=2 . .nops (lpol) )] ) , [dtvars] ) 

end do: 

# The only term left is the correct one. 
ode : =op(lpol) : 

# Convert to recurrence. 
REC : =dif f eqtorec ( 

{applyopr(ode, F(t||k), At), F(0)=1}, F(t||k), a(n)): 

# Calculate some terms. 

GRAPH : =rectoproc (REC , a(n) , list) (20) : 
[seq(GRAPH(10) [i]*(i-l) ! ,i=l. .20)] ; 



[1 , , , , 1 , , 70 , , 19355 , , 11180820 , , 1 1555272575 , , 

19506631814670 , , 50262958713792825 , , 187747837889699887800 , 0] 
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